Disorder effect on 3-dimensional Z 2 quantum spin Hall systems 
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In this paper, we address ourselves to the nonmagnetic disorder effects onto the quantum critical 
point which intervenes the 3-dimensional Z2 quantum spin Hall insulator (topological insulator) 
and an ordinary insulator. The minimal model describing this type of the quantum critical point is 
the single-copy of the 3+1 Dirac fermion, whose topological mass m induces the phase transition 
between the topological insulator and an ordinary one. We first derive the phase diagram spanned 
by the mass-term m, chemical potential /i and strength of the disorder within the self-consistent 
Born approximation. By way of this, we find a finite density of state appears even at the zero- 
energy and at the phase transition point, i.e. m = /i = 0, if the strength of the disorder potential 
exceeds some critical value. To infer the structure of the low-energy effective theory around these 
zero-energy states, we further calculated the weak localization (WL) correction to the conductivity. 
To be more specific, we have found that the diffuson is dominated by the charge diffusion mode and 
parity diffusion mode. While the charge diffusion mode always carries the diffusion pole, the parity 
diffusion mode becomes massless only at m = 0, but suffers from the infrared cutoff for non-zero m. 
Corresponding to this feature of the diffuson, the Cooperon is also composed of two quasi-degenerate 
contributions. We found that these two give rise to the same magnitude of the anti-weak-localization 
(AWL) correction with each other at m — 0. As a result, when the topological mass m is fine-tuned 
to be zero (but for generic /1), the AWL correction becomes doubled (quantum correction doubling). 
Based on this observation, we will discuss the possible microscopic picture of the "levitation and 
pair annihilation" phenomena, recently discovered by Onoda et ai2&. 
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I. INTRODUCTION 

Physics of spin transport has been a matter of inten- 
sive research in recent years. One of the topics of current 
interest is the spin Hall effect. This has been originally 
proposed theoretically^ and later followed by various 
experimental results 3 ' 4 -. The research on the spin Hall 
effect opens a new field of Hall effects in time-reversal in- 
variant systems. This has also led us to a new concept of 
the quantum spin Hall effect, which is the natural "spin" 
extension of the quantum Hall effect^^. In the quan- 
tum spin Hall effect in two dimensions, the bulk is gapped 
while there are gapless edge states carrying spin current. 
In this case, the external magnetic field is zero, while the 
spin-orbit coupling acts as a "spin-dependent magnetic 
field", giving rise to the effect analogous to the quan- 
tum Hall effect. Such insulators showing the quantum 
spin Hall effect are characterized by the Z2 topological 
number—. We shall call them as Z2 quantum spin Hall 
insulators (QSHI). 

The simplest system for the two-dimensional (2-d) Z2 
QSHI is realized as a superposition of the wavefunctions 
of two quantum Hall subsystems ii&I, one with spin up 
and the other with spin down, having opposite Chern 
numbers. The system respects not only the time-reversal 
(T) invariance but also the spin-conservation. Such an 
insulator supports the same numbers of right-moving up- 
spin edge states and left-moving down-spin edges. 

This Kramers pair of chiral edge states is often called 
as the helical edge state. By its construction, the num- 



ber of this Kramers pairs of edge states correspond to 
the Chern integer associated with its bulk wavefunctiorA 
The T symmetry guarantees the double degeneracy be- 
tween right-moving up-spin and left moving down-spin 
states. Thus, the stability of each Kramers pair is sup- 
ported by this T-symmetry. 

Spin non-conserving (but T invariant) perturbations, 
however, introduce level repulsions between two differ- 
ent Kramers pairs. Namely, they usually let two pairs 
annihilate with each other, and open a gap. Accord- 
ingly, in the presence of generic spin-non-conserving per- 
turbations, those wavefunctions having even numbers of 
Kramers pairs, reduce to trivial insulators, which have 
no gapless edge states ^ 10 ' 11 ? 12 . Meanwhile, wavefunc- 
tions having odd numbers of pairs still can have one ac- 
tive helical edge mode. The latter is dubbed as the Z 2 
quantum spin Hall (topological) insulator. Thus, stabil- 
ity of such a gapless edge state is protected only by the 
T-symmetry, while does not require spin-conservations 
anymore^ 1 ^ 1 ^ 3 .^. 

The three-dimensional (3-d) versio n 15 i 16 i 17 ' 18 i 19 i 21 of 
the Z2 QSHI carries same characters as that of 2-d does. 
The 3-d Z 2 QSHI also allows any spin-nonconservering 
perturbations, while always requires the T symmetry. Si- 
multaneously, however, it is not a mere extension of the 
2-d Z 2 QSHI, in a sense that 3-d Z 2 QSHI has no U(l)- 
analogue of QSHI. Namely, they support a 2 + 1 massless 
Dirac fermion as its surface state 15 i 19 , instead of a helical 
edge state. In the 2-d surface Brillouin zone, say k x -k y 
plane, this massless Dirac fermion has a spin which de- 
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FIG. 1: A schematic phase diagram for the quantum critical 
point intervening the Z2 QSHI and an ordinary insulator, m 
is a system parameter driving the phase transition. When 
m = 0, the system is in a critical phase. In the Fu-Kane-Mele 
model, it corresponds to the relative strength of one out of 
the four NN transfer integrals emitting from a single site. 



pends on the (surface) crystal momentum. It is clear that 
such an insulator cannot be adiabatically connected into 
a composite of two spinless wavefunctions. In such Zi 
QSHI, the T symmetry therefore guarantees the mass- 
less nature of each 2 + 1 surface Dirac fermion. 

Zi QSHI always has a quantum critical point at its 
phase boundary to any ordinary insulators, in both 2-d 
and 3-d. For example, from the 3-d tight-binding model 
proposed by Fu, Kane and Melei^, we can explicitly see 
this; when a certain T symmetric parameter is varied 
in their model, 3d Zi QSHI is driven into an ordinary 
insulator, latter of which does not support any surface 
states (see Fig. [T]). Observing this, a following question 
naturally arises; during this tuning, the T symmetry is 
always preserved, so that the massless nature of the sur- 
face Dirac fermion is supposed to be protected by this. 
At the same time, however, this 2 + 1 surface fermion 
should have become "massive" , when a system enter an 
ordinary insulator phase. Thus, one might ask how this 
single surface Dirac fermion could acquire a finite mass, 
with keeping the T symmetry? 

The answer is simple; we have two sample boundaries, 
say z — +L and z = —L. Each boundary supports one 
2+1 surface massless Dirac fermion respectively. They 
are localized at each boundary, when the bulk gap is 
sufficiently large. In such a situation, a mixing between 
these two 2 + 1 surface massless Dirac fermions is tiny, i.e. 
0{e~ L /^) with £ being the localization length. However, 
when a system becomes close to the quantum critical 
point, a mixing between these two surface states becomes 
substantial, with increasing £. When a bulk eventually 
reaches the quantum critical point, two surface massless 
Dirac fermions readily communicate via extended bulk 
states. Thus, they generally annihilate in pairs, just as 
in those insulators having even number of 2 + 1 surface 
Dirac fermions at one boundary. 

This simple picture in the clean limit raises the follow- 
ing non-trivial speculations about the disorder effects on 
the Zi QSHI. Suppose that T-symmetric random poten- 
tials are introduced in the topological insulator phase. 
When the corresponding bulk gap is sufficiently large, 
we could begin with two separate bands. The scaling ar- 
gument in 3-d22 tells us that each band should always 
have two mobility edges, respectively (see Fig. HJb)). 
Namely, there is no delocalized bulk-wavefunction near 
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FIG. 2: A schematic picture of the density of state and mobil- 
ity edges, where hatched region corresponds to the extended 
state. Inset represents the energy dispersion as a function of 
surface crystal momentum, where the red line corresponds to 
the surface state at z — ±L. (a) m — 0; at quantum critical 
point, (b) m > 0; in the topological insulator phase. 



the zero energy. Starting from this phase, let us change 
some T-invariant model-parameters, so that a bulk tran- 
sits from this topological insulator to an ordinary one. 
From the argument in the clean limit, one can then ex- 
pect that a delocalized bulk-wavefunction should emerge 
at the zero- energy region at the quantum critical point, 
i.e. {i — m — (see Fig. dja)). If it were not, the 
two surface states localized at the two sample boundaries 
could not communicate at all and they could not annihi- 
late with each other. As a result, the system was unable 
to smoothly enter an ordinary band insulator, since the 
latter one does not support any surface state at all. 

To put this reversely, the existence of the quantum 
critical point (QCP) having extended bulk wavefunctions 
is always required, whenever this critical point separates 
an ordinary insulator and the topological insulator. This 
is because these two insulating phases support different 
numbers of Kramers pairs of surface states. Moreover, 
provided that these surface states are stable by itself, this 
QCP should be also stable, however small the density of 
state (DOS) at the zero energy is and however strong the 
disorder strength is. Otherwise, the topological insulator 
could be adiabatically connected into an ordinary band 
insulator, which contradicts the different Zi topological 
numbers for the two phases. 

In this paper, we will uncover several novel features as- 
sociated with the non-magnetic disorder effects onto this 
topological quantum critical point. The organization of 
this paper is summarized as follows. In the next section, 
we will briefly review the effective continuum model for 
the quantum critical point intervening the Zi topological 
insulator and an ordinary insulator. The effective model 
is known to be described by the 3+1 Dirac fermion, whose 
mass term brings about the topological quantum phase 
transition. Namely, when the mass term changed from 
positive to negative, a system transits from the topolog- 
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ical insulator to an ordinary insulator. As such, we call 
this mass term especially as the topological mass term. 
Based on this effective model, wc will next introduce var- 
ious types of the on-site random potentials respecting the 
T-symmetry. Note that, in this paper, we restrict our- 
selves to T-symmetric cases and exclude magnetic impu- 
rities, because, in the absence of the T symmetry, the 
two phases are no longer topologically distinct. 

Based on the self-consistent Born approximation, we 
first work over the single-particle Green function in the 
section III. The phase diagram spanned by the (bare) 
chemical potential [i, (bare) mass term m and strength 
of the disorder a is derived. In particular, at the critical 
point, i.e. m = 0, we found some critical value of the 
disorder strength, a Cl above which the zero-energy state, 
i.e. [i = 0, acquires a finite life-time t; 
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Since the density of state in our model is always pro- 
portional to the inverse of the life-time (see below) , non 
zero r" 1 simply means that a system is in a compressible 
phase. 

When a finite but small topological mass m is intro- 
duced for a > a c , the life-time t and the renormalizcd 
mass to becomes as follows; 
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where tq is given as a function only of a via Eq. (fTj). 
Thus, when the bare topological mass exceeds the critical 
value m c = 2t ~ 1 , the density of state vanishes, so that 
a system enters an incompressible phase. This gapped 
phase can be adiabatically connected into band insula- 
tor phases in the clean limit. Accordingly, we will reach 
the phase diagram for the zero-energy state as depicted 
in Fig. [3] In the section III, wc also describe the be- 
havior of the one-particle Green function for a finite /i 
(see below), in which we observe that the compressible 
phase (not necessarily metallic phase) always intervenes 
the topological insulator phase and an ordinary insulator 
phase as in Fig. [3] 

Focusing on this intervening compressible phase, espe- 
cially for a < a c , we will derive in the section. IV the 
diffuson, Cooperon and the weak localization correction 
to the electric conductivity. We will first observe that the 
diffuson is composed of two quasi-degenerate low-energy 
modes; 
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(see Fig. (jlll) or eq. ([74]) for the definition of T d (q,Lu)). 
The first term participates in usual charge diffuson mode, 
and therefore always has the diffusion pole, i.e. [u> + 
iDq 2 ]^ 1 (lu and q stand for the frequency and momentum 




ordinary insulator topological insulator 



m 



FIG. 3: A schematic phase diagram for /i = 0. A blue shaded 
region corresponds to a compressible phase, which separates 
two gapped phases, i.e. an ordinary insulator and the topo- 
logical insulator. This phase boundary for a > a c is given by 
1 _ 2s = S ArcTan [-2-1. 

a 2 L tn J 



of the density fluctuation respectively). The other low- 
energy mode, however, becomes massless only in the ab- 
sence of the topological mass m. Namely, its low-energy 
and long wavelength behavior is generally truncated by 
the infrared cutoff T t ~p OJ while this infrared cutoff reduces 
to zero at to = (but generic /i), i.e. T t ~p cx to 2 . 

Physically speaking, this second mode describes the 
diffusion of the parity density degree of freedom, which 
becomes a conserved quantity of our effective hamilto- 
nian at to = 0. Namely, the parity-density correlation 
function exhibits the diffusion pole structure at the crit- 
ical point (to = 0), while it becomes massive in the pres- 
ence of the finite topological mass. Consequently, the 
diffuson acquires one additional low-energy, i.e. the 2nd 
term of eq. ([3]) , into which the information of this parity- 
density correlation function is separately encoded. 

When the hole-line of the diffuson time-reversed, these 
two-mode features are transcribed into the Cooperon: 
the Cooperon thus obtained is also composed of two 
quasi-degenerate dominant contributions; 
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(see Fig.fTTIand its caption for the definition of U coop (k + 
k',cu)). In section. IV, we will see that, at the criti- 
cal point (to = 0), these two contributions give rise to 
the same amplitude of the anti-weak-localization (AWL) 
correction to the electric conductivity. When the finite 
topological mass is introduced, however, the Cooperon 
associated with the parity mode channel becomes ineffec- 
tive, since its backward scattering behaviour is truncated 
by the cutoff r t ~p . Meanwhile, the Cooperon obtained 
from the charge mode channel remains effective, even in 
the presence of finite to. As a result, the AWL correction 
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at the critical point becomes precisely halved, on intro- 
ducing the finite topological mass (quantum correction 
doubling) . 

In terms of this novel behaviour of the parity diffu- 
sion mode and that of the corresponding AWL correction, 
we will argue in the section. V the possible microscopic 
mechanism of how the delocalized bulk-wavefunction 
emerges at the critical point, i.e. m = fJ, = 0. To be more 
specific, we expect that the parity diffusion mode men- 
tioned above generally becomes massless, when a system 
transits from the topological insulator side to the ordi- 
nary insulator side. Assuming that this is the case, we 
will attribute the emergence of the extended bulk wave- 
function to the AWL correction obtained from this parity 
mode channel, i.e. the second term of eq. ((H). For the 
systematic understanding, however, one generally needs 
to go beyond the theoretical approach employed in this 
paper. Several open issues will be also discussed in the 
section V. 

A number of appendices describe other topics useful in 
understanding the main text in more detail. For clarity of 
the explanation, we have presented the results only in the 
case of chemical-potential type disorder in the text. The 
study in the presence of general T symmetric disorders 
becomes more cumbersome. But the basic feature such 
as the phase diagram is expected to be same. In the 
appendix A, we will describe how the one-particle Green 
function at the zero-energy state behaves in the presence 
of these general T-invariant random potentials. 

Our weak-localization calculation described in the sec- 
tion. IV is the controlled analysis, when it comes to the 
weakly disordered region, a < a c . Namely, for this 
parameter region, one can confirm self-consistently the 
coupling constant 1/kpl = 1//j>t to be sufficiently small 
around fi = (see eq. ((70)) ). For a > a c , however, this 
coupling constant generally diverges toward fi = 0, only 
to make the weak localization calculation (and scB cal- 
culation) an uncontrolled analysis. Thus, as the com- 
plementary analysis for this strongly disordered region, 
a > a c , we employed the mode-mode coupling analysis 
in the appendices. B-D. By taking into account the quan- 
tum interference effect due to the Cooperon terms, this 
theoretical framework gives us the gap equation for the 
diffusion constant. Main results in section. IV such as 
the quantum correction doubling are also supported by 
the analysis in the appendices. B-D. 



II. EFFECTIVE CONTINUUM MODEL AND 
DISORDER 

A. Effective continuum model 

We consider a system with both T-symmetry and the 
spatial inversion (Z)-symmetry. Under this symmetry re- 
quirement, Murakami et al. recently derived the minimal 
model for an arbitrary quantum critical point interven- 
ing the topological insulator and an ordinary insulator 



on a quite general groun d ' 25 . It turns out to be always 
described by the 3+1 Dirac fermion given as follows; 

r 3 

Ho = d 3 r^(r){j2%{-^)-m%}^(r), (5) 

where m corresponds to the topological mass term. 
Without loss of generality, one can regard the topological 
insulator phase to be ra > and an ordinary insulator 
phase to be m < (see Fig.[T]). To see that "m" actually 
endows this Dirac fermion with a mass, we note that fol- 
lowing five 4x4 7-matrices are anticommuting with one 
another; 

7i = fry (g) 1, 72 = fr z <8> s x , 73 = fr z (g) Sy, 

74 = frz ® S z ,7 5 = fr x <g> 1. 

The matrices fr^ and are Pauli matrices, representing 
the (generalized) sublattice index, and the spin index, re- 
spectively. In terms of these Pauli matrices, we will take 
the T operator as is y K with K being the complex con- 
jugation. Meanwhile, the X operator will be taken as fr x . 
It follows from these conventions that 71,2,3,4 are T odd 
and X odd, while 75 is T even and X even (see table. I). 
Together with the property that — id^ is T odd and X 
odd, we can easily see that our Hamiltonian in the clean 
case is indeed T even and X even. This guarantees the 
Kramers degeneracy at each fc-point, irrespectively of the 
topological mass m. We also note that eq. (f5]) is indeed 
the low-energy effective continuum Hamiltonian for vari- 
ous lattice model recently discussed in literature o 15 ! 19 ' 20 . 
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72 = fr z ® S x 






73 = fr z ® Sy 






74 = fr z (g> S Z 






75 = fr X <g> 1 
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TABLE I: Dirac operators and their symmetries. 

Generally speaking, we can enumerate all Hermite ma- 
trices possible in this spin-sublattice space. Namely, us- 
ing the commutator between these five Dirac matrices, 
we have other 10 = 5C2 associated Dirac matrices; 

% = ^[7i,7j] = -HiTr (6) 

We can further classify these 10 matrices into two classes; 
one is T invariant (even) matrices and the other is T odd. 
Since the five Dirac matrices are always even under X ■ T, 
these 10 associated Dirac matrices are by construction 
odd under X-T. Thus the symmetries of these 10 matrices 
can be summarized as in Table II. 
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TABLE II: Dirac associated operators and their symmetries 



let us introduce T-symmetric "on-site type" random 
potentials as generally as possible; 

r 4 

Wimp = / dr^(r)^v jQ + vr o % + ^2v j5 j 35 ^tp(r), (7) 

where all the 6 components of the vector v(r) are real- 
valued functions of r. Then, each single-particle eigen- 
state of Ho + Himp always has a Kramers pair state; 



(4>(r)\ = i® (-*)$„!#»•)>. 



(8) 



Namely, \(f>(r)) and \4>{r)) are degenerate and orthogonal 
to each other. Noting this, one can see that the retarded 
(advanced) Green function observes the following relation 
in each ensemble; 



G«W(r,r';/i) ee £ 



I0n(r))(0„(r')| 
fj, - e„ ± i<5 



(9) 



B. Disorder averages, spatial inversion symmetry, 
rotational symmetry and engineering dimension 

As usual, we will take the quenched-average of these 
T invariant impurities at the gaussian level; 



P[v] = j i d 3 rd 3 r' 

j',to£{0,5,15,--,45} 

[^^(rjy ,m)Vj(r)v m (r'), 



(10) 



with a proper normalization factor J\f and real-valued 
symmetric matrix A. For simplicity, an "on-site type" 
correlation will be assumed; 



A{r,j\r',m) = A jm S 3 {r-r'). 



We also suppose that the translation symmetry and 
the spatial inversion symmetry are recovered after these 
quenched averages; 

G R{A) (r,r';[i) = G R{ - A \r + b,r' + (13) 
qR(A) (r> r /. m) = a x ®\- G R W {-r, -r'; M ) • & x ®(l4) 



Then, the latter symmetry, i.e. eq. (| 14(1 . prohibits any 
matrix elements between 70,5 and 7^5 (j = 1, • • • , 4) in the 
right hand side of eq. (fl2|) . Namely, the 6x6 matrix A 
in its right hand side takes the following form; 



(15) 
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A 50 
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A a 
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U515 



14515 



11545 



11515 



(16) 



This is because 70,5 are even under 2, while 7^5 (j = 
1,2,3,4) are odd. In order that the gaussian integral in 
eq. (|10p converges, all the eigenvalues of A have to be 
positive. Accordingly, the matrix elements described in 
eqs. (I15H16P must obey the following inequalities; 



A00A55 > A05A50 — A 05 , 
Aoo + A 55 > 0, TrA a > 0, 



(17) 



We can study the effects of these general T-invariant 
"on-site type" random potentials, without any further as- 
sumptions. As will be partly shown in the Appendix A, 
however, such an analysis becomes very cumbersome and 
lengthy. Thus, we henceforth consider only the chemical 
potential type disorder Aoo > because it is expected to be 
dominant among various types of disorder. Those who 
are interested in the effects of other components such as 
Ao5,A55 and A a may consult the appendix A. In Ap- 
pendix A we have studied the effect of the T-reversal 
invariant "on-site type" disorder on a general ground, 
focusing on the zero-energy wavefunction at the critical 
point. 

Being translationally invariant as in eq. (|13p , the aver- 
aged Green functions can be readily fourier-transformcd 
by the use of the crystal momentum k. The resulting 
Green functions can be expanded in terms of Dirac ma- 
trices and its associates; 



G R (k,fi) 



E 



(18) 



je{0,l,-,5,15, — ,42} 



(11) c 



Fj(/c, /it) stands for some complex- valued function of A; and 
/.t. In this momentum representation, T and X invariance, 
i.e. eq. © and eq. (TL4"|) , read as follows; 



(12) 



a x ®i-G R ( A \k,»)-& x 
i®s v -G R ( A \k,(i) • 1 



i = G ^)(-fc, M ), (19) 
h = {G fl( - 4 )} t (-fc, M )(20) 
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These two symmetries require that Fj = i ) ... ) 4(fc, /z) are odd 
functions of k, Fo,5 are even functions of k, and also that 
Fjj = for i =/= j and i, j — 1, • • • , 5. Namely, the retarded 
and advanced Green functions are given only in terms of 
the anti-commuting Dirac matrices; 



G R (k^) = F (fc,/i) l + ^F^/x) 



7^, 



7a»- 



(21) 



(22) 



In addition to the T-symmetry and X-symmetry, the 
pseudo-spin rotational symmetry is also recovered after 
the quenched average. This is because only the chemical- 
potential type disorder Aoo is considered now. Specifi- 
cally, the 1-point Green function after the quenched av- 
erage respects the simultaneous rotations of the spatial 
coordinate and the pseudo-spin coordinate; 

C> n ,0-G K(A) (fc^)-t>L = G R(A) (Rn,<t>k,V>), (23) 



n,tp 



/i, v and p above run over 1, 2 and 3. R n .<f, in the right 
hand side stands for the spatial rotation around the vec- 
tor n by the angle (j>. When combined with eqs. (|19l20p . 
this rotational symmetry further restricts the form of 
the Green functions. For example, the coefficient of 74 
should be an odd function of k due to eq. (fl"9|) . while it 
should be an even function of k because of eq. (|2"5|) . As 
such, Green functions cannot contain 74-component, un- 
der these two symmetry requirements. Moreover, eq. (|23[) 
by itself compels Fi,2,a(k, (i) to be transformed as a vec- 
tor under the rotation in the fc-space; 



c 3 fc 2 fc M 



(24) 



with fi = 1, 2, 3. 

So far, we have imposed several generic symmetries 
such as T-symmetry and X-symmetry on the Green func- 
tion after the quenched averaged. As a result of this, the 
Green function is given only in terms of the Dirac matri- 
ces. Since these 5 Dirac matrices are all anticommuting 
with one another, the inverse of the Green function can 
be easily calculated, 

5 

G^- l {k,ti)= F (fc, M ) l + J2 F »( k >») >, (25) 



_Fo = F 



F = L 4 , F v = —, • (26) 

r 2^=1 r M 



Correspondingly, the inverse of the bare Green function 
is given as follows; 



Gf~\k^) = (ji + iS)i- J2 k *^ + 
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A=l,2,3 



A=0,---,5 



(27) 



Based on these simplifications, we will derive in the 
next two sections the electronic property of the disor- 
dered single copy of 3 + 1 Dirac fermion described by 
eq. (|5|). Before finalizing this section, however, it would 
be appropriate to summarize the engineering dimension 
of the various quantities introduced in this section. Com- 
paring the impurity hamiltonian with the pure hamilto- 
nian, one can first see that 



m,fe At ,/i,f /1 ,F M ,u i ~ [L l ], Fn~[L], 



(28) 



where L denotes the dimension of a length. Out of this, 
we can further figure out the engineering dimension of 

A im , 



(29) 



by requiring P[v] in eq. (JTTJ) to be dimensionless. 



III. SELF-CONSISTENT BORN 
APPROXIMATION 

The self-consistent Born (scB) approximation simply 
equates the right hand sides of the following two; 



S"(fc, M ) = G R '- 1 -G^ 
5 

= (f - F ) 1 + £)(f„ - F„) %, 



VR.-1 



(30) 



E"(fc,/x)=A o / d i k'G tt {k\ l x) 



= A 00 / d 3 fc'{F (fc', M )7o + F 5 (fc', M )7 5 }(.31) 

J0<\k\<A 

We have already omitted those terms proportional to 
Fi,2,3,4(fc', fj) m the integrand of eq. (f3"l"j) . since they are 
odd functions of k'. Comparing the coefficients of each 7 
matrix in eq. (|5D|) and eq. (|3"Tj) . we can make the closed 
coupled equation for Fq and F5; 



d A k 



/ 0< I fc I < A 
-Ann / d 3 k 



= fo - F , (32) 



0<|fc|<A 



F 2 - F 2 - k 2 



(33) 



by the use of eq. (|26|) . We have already used the following 
relations also: 



Fl,2,3 = 'l,2,3 — — ^1,2,3) F 4 = f4 = 0. 



(34) 



These integral equations in eqs. (|32j|33p clearly depend 
on the ultraviolet cut-off A. Thus, rescaling the momen- 
tum by this cut-off A, let us introduce the dimensionless 
quantities, instead of F M , f M , and Aoo- Eqs. (|28ti29|) indi- 
cate that they should be rescaled in the following way; 



Fo,5 - * Fo,5 = Fo.sA , 

fo,5 — ► /o,5 = fo,5A , 

A 00 —> a = 27rA 00 A. 



(35) 
(36) 
(37) 
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The factor 2-k in the definition of a is just for later con- 
venience. In terms of these dimensionless quantities, the 
above coupled non-linear equations become; 

(l + aG)F = f =n±i5, (38) 
(l-aG)Fg = / 5 = 77i, (39) 

where fj, and m in the right hand side are supposed be 
also normalized by A -1 . G used in the left hand side was 
also made dimensionless; 

"■'Lfiib^ (40) 

(a + ibfsFS-F?. (41) 
Eqs. (|38H4ip thus determine F Q and F 5 as a function of 



their bare values; /o and /g. i 7 ^ thus obtained should be 
by definition much smaller than the "ultraviolet cut-off" 

i; 

F„ < 1. (42) 

This also leads to a, b <C 1. In the fallowings, we will 
frequently take full advantage of their smallness, which 
is always self-consistently verified later (see below) . 

In the next subsection, we will present the solution 
of this coupled equation for general fi and m. Before 
doing this, however, it would be appropriate to express 
the imaginary part and real part of G in terms of a and 
6, so that one can roughly estimate these two quantities 
in small a and b. The real part and the imaginary part 
of G read as follows; 



ReG 
ImG 



log 



\l-a) 2 + b 2 



.(1 



b 2 



-2 log 



b 2 



fr^ArcTan 
a ( ArcTan ' 



ArcTan 

•n 



ArcTan 



(43) 
(44) 



Observing these two, please notice that the final two 
terms in eq. (|43M4j) . which are proportional to ArcTan, 
are nothing but the pole contribution. Namely, the limit 
b — > reduces them a finite constant with its sign iden- 
tical to that of b, e.g. 

ArcTan^ 1 - + ArcTan^ 1 ~l ° ) — > 7rsgn(&), 

where one should also note that a,b <C 1. Bearing these 
in mind, one can then evaluate the leading order of ReG 
and ImG with respect to small a and 6; 

ReG = -2 + tt|6| + 0(a 2 ), (45) 
ImG = -sgn{b)Tra + 0(ab). (46) 

Namely, the second member of eq. (I45|) and the first mem- 
ber of eq. (|46|) are nothing but the pole contributions 
mentioned above. 

A. Solution for scB equations 

1. m — n = case 

For the warming-up, consider first the case with /i = 
m = 0, i.e. the zero-energy state at the critical point. 
Equations (138H4ip have three types of solutions; 

(i) : ^ = ^ = 0, (47) 

(ii) : 1 + aG = n F 5 = 0, (48) 

(iii) : 1 - aG = n F = 0. (49) 



Observing the estimates given in eqs. (I45ti46|) , please no- 
tice that type- (iii) solution cannot be satisfied for a > 
and a,b <^ 1. Thus, we will ignore this henceforth. 

The type-(i) solution is always trivially satisfied. This 
solution indicates that the zero-energy state is not renor- 
malized at all by the disorder, Fq = fo = 0, F$ = /s = 0. 
Thus, it describes the diffusionless zero-energy state. 

The type-(ii) solution is a non-trivial solution, which 
turns out to describe the diffusive zero-energy state. To 
see this, let us begin with the first condition of eq. (|48[) . 
i.e. I + aG = 0. The imaginary part of this gives ImG = 
0, which is satisfied either when a = or when 6 = and 
\a\ > 1 (see Fig. [4} . Since a, b <C 1 as noted earlier, the 
only physical solution satisfying ImG = is thus a = 0. 
The remaining condition, 1 + aReG = 0, becomes then 
simple; 

6ArcTan[6" 1 ] = 1 - (50) 

Since F§ — gives a + ib = F , a and b thus obtained 
stand for the renormalized chemical potential Ti and the 
inverse of the lifetime r _1 , respectively. Accordingly, 
the type-(ii) solution simply denotes that the zero-energy 
state acquires a finite lifetime r, while its chemical po- 
tential is free from renormalizations; 

11=0, r^ArcTanM = 1 - — . (51) 
L J a 

Namely, for a > a c = h, r _1 can take a finite value. 

For a weak disorder region (a < a c ), Eq. (|51|) cannot 
be satisfied for any r. Thus, the only solution therein is 
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(a) 



+ 


1 


-I 


+ 



a 




ImG 




ImG 



FIG. 4: ImG as a function of a and 6. (a) The sign of ImG, 
which is an odd function both in a and in b. "+(—)" stands 
for the sign of ImG at the 4 regions, i.e. a, b > 0, a > > b, 
b > > a and > a, b. The bold line which runs from ( — 1, 0) 
to (1, 0) denotes a sort of the branch cut. Namely, ImG jumps 
from -2n 2 a to +2iv 2 a (from b = +0 to b = -0). (b) A side 
view plot of ImG only for b > 0. (c) ImG = is satisfied 
either when a — 0, or when 6 = and \a\ > 1. 



the type-(i) trivial solution. On the other hand, both the 
type-(i) solution and type-(ii) solution become possible, 
above this critical disorder strength (a > a c ). In the 
next three paragraphs, we will determine which solution 
is physically sensible for a > a c . 

To do this, we will extend these two solutions into a 
small but finite /i region. Namely, by seeing how this 
chemical potential will be renormalized for each case, we 
will judge which solution is the physical one for a > a c . 
Recall first that F5 = in either case. Thus, a and 
b correspond to ~p and r _1 respectively, so that a and 
b should be an odd and an even function of the bare 
chemical potential fj, respectively. 

Bearing these in mind, let us extend the type-(i) so- 
lution into a small fx- region first. Namely, keeping the 
leading order in small fx, we can evaluate the real part of 
eq. J5BJ); 



(1 - 2a)a + 0{fi 3 ) = fx, 



(52) 



where we used a cx 0{fx) and b cx <D{fi 2 ). Thus, the 
renormalized chemical potential is estimated up to C(/i) 
as follows; 



fi = a 



1 - 2a 



(53) 



while b will be determined up to (0)(/i e ) from the imag- 
inary part of eq. ([38]) : 



Q7T 



(1 -2a) 3 



fx 2 + 0(» 4 ). 



(54) 



This solution indicates that the negative eigen-energy 
state and the positive eigen-energy state are inverted en- 



ergetically for a > a c ; sign/I = — sign//. This is, however, 
clearly unphysical at least for small fi. 

When the type-(ii) solution is extended into a small fi- 
rcgion, a similar algebra gives us the following expression 
for the real part of eq. (f3"5)) up to O(n); 



(1 - 2a)a + 4ar 1 ArcTan[r]( 



(55) 



In this equation, we have already made implicit those 
contributions proportional to C(/i 3 ) and C(r _2 /z) while 
keeping those proportional to C(r -1 ^) explicit. Please 
also note that we have used a cx 0{n) and b = r _1 + 
C(/i 2 ). Namely, contrary to the type-(i) solution, b con- 
verges to a non-zero r _1 at the leading order in small 
/i. As a result of this, eq. ([55]) has acquired an addi- 
tional C(/i)-term, i.e. 4ar _1 ArcTan -1 [r] • a, which was 
absent in eq. (|52[) . This additional term makes the sign 
of JI to be same as that of [i. Namely, by the use of 
l-2a = -2UT- 1 ArcTan" 1 [t] +0(fi 2 ), eq. (55]) leads us 
to; 



(2a- l)a + C(^ 3 ) = fi. 



(56) 



Out of this, one can evaluate the renormalized chemical 
potential up to O(fi) as follows; 



fj, = a 



fi 



2a - 1 



(57) 



whose sign is clearly same as that of the bare one for 
a > a c ; sgn/Z = sgn/i. Observing these two distinct 
behaviors for the finite fi region, i.e. eq. ([55]) and eq. ([57]) . 
we conclude that, for a > a c , the type-(ii) solution is the 
physically sensible solution, while the type-(i) solution is 
an unphysical one. 

To summarize so far, the physical solutions obtained 
at m = fj, = are two-fold, depending on the disorder 
strength a. When a < a c = 1/2, the type-(i) trivial 
solution is the only possible solution; 



(i) : F — F 5 = for a < a c 



(58) 



This means that the electronic state at the zero-energy 
is free from the disorder up to a certain critical disorder 
strength. 

On the other hand, when its strength exceeds this crit- 
ical value, i.e. a > a c , the type-(ii) solution should be 
adopted; 

(ii) : F = it' 1 , F 5 = for a > a c . (59) 

This solution means that the electronic state at the zero 
energy acquires a finite lifetime r defined by eq. (15ip . 



2. n = and finite m case 

Let us introduce a finite topological mass m into 
eq. ([55]) and eq. ([B"9"]) respectively, with the chemical po- 
tential fi being still zero. We will first argue that the 
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solution of eqs. (|38H41|) in the presence of the finite mass 
is uniquely determined for a < a c . Such a solution reads; 

F = 0, F 5 =m, (60) 

where to is given as a function of the bare mass; 

m{l + 2a - 2amArcTan[m~ 1 ] } = to. (61) 

A typical behavior of to as a function of the bare mass is 
depicted in Fig. 02b). 

To see that eqs. (|60M6ip is the only possible solution 
for a < a c , let us begin with the real part of 1 + aG 
appearing in eq. (138p . In the case of a < a c , it is always 
positive definite for any a < 1. As such, we must take Fq 
to be zero, to satisfy eq. ([3"5| . This leads to F& = b — ia. 
Using this, consider next the imaginary part of eq. ([39]); 

(1 - aReG) ■ a + almG -6 = 0. (62) 

Observing the leading-order estimates of RcG and ImG, 
i.e. eqs. (|45M6j) . one can further see that eq. (|62|) 
uniquely leads to a = 0. The remaining condition, i.e. 
the real part of eq. ([39]) , then becomes simple; 

(1 - aReG) • b = (l + 2a - 2a6ArcTan[V 1 ] ) ■ b = m. 

Now that (Fo,F$) = (0,6), this equation is nothing but 
eq. (|61[) . when b replaced by to. 

Let us next consider the case of a > a c . The solution 
of eqs. (|38H4ip in this case is two fold; we have a certain 
critical mass value m c , which is given as a function of a; 

m c = 2T- 1 , r^ArcTanM = 1- — . (63) 

a 

When the topological mass is less than this critical value 
(to < to c ), the solution of eqs. (I38M1|) becomes; 

F Q = +Wt- 2 -m 2 , F b = m = —, (64) 

where r was already defined in eq. (|63p. On the other 
hand, when the topological mass exceeds this critical 
value (to > to c ), the solution becomes eqs. (I60H61P again. 

To see that eq. fl64} is the solution of eqs. (|38I41|) 
for a > a c and to < to c , take 1 -I- aG = first, so 
that eq. (f3"5| is satisfied. By the use of the same argu- 
ments described above eq. (|50p , this immediately gives us 
(a,b) = (0,t _1 ), with r being defined by eq. ([ST]) . Since 
1 - aG = 2, eq. ((3SJ) leads to F 5 = to/2. Thus, using 
these two things, we obtain Fq out of eq. (|4"Tj) . which 
is nothing but eq. (|64[) . When to exceeds m c = 2r _1 , 
eq. (|64p becomes an unphysical solution in a similar way 
as the type-(i) solution in the previous subsection did for 
a > a c ; 

I n 77 TO 

Fq = ± V TO — T , i*5 = TO = — . 

Instead of this, it turns out that we should adopt the 
other solution for to > to c , i.e. eqs. (|60H6ip . 




0.0 0.2 m 0.4 




FIG. 5: (a); A schematic phase diagram at /1 = 0. The 
white region corresponds to the incompressible phase, where 
no finite DOS exists at /1 = 0. The blue region corresponds to 
the compressible states, where a finite DOS exists at /1 = 0, 
i.e. eq. (|64p . (b) The renormalized mass m as a function of 
the bare mass m for a — 0.25 < a c . (c) m as a function of 
m for q = 0.75 > a c . There exists the critical value of the 
bare mass m, below which m = m/2, and above which m is 
determined by eq. (I6ip . These two values coincide with each 
other precisely at m = m c . 



A typical behavior of to in the case of a > a c is de- 
picted in Fig. [51(c), where these two solutions, i.e. eq. |64|) 
and eqs. (|60H6ip . are indeed connected continuously at 
to = to c . Since eq. (|64p always supports a finite den- 
sity of state (DOS), we can regard that the compressible 
phase extends over a > a c and to < m c . On the other 
hand, eqs. (|60H6ip do not support any finite DOS. As 
such, we can regard that an incompressible phase extends 
over a < a c or to > to c (Fig. [Sfa)). 
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H>0, m<0 ^ \i>0, m>0 



1.5X16 



(Fo,FJ,-F5,-F5 ) : retarded 
(Fq -Fq ,-F 5 ',F5 ) : advanced 



(-Fo',Fo,-F 5 ',F5) : retarded 
(-Fo — F ",-F5,— F5 ) :advanced 



|i<0, m<0 



(Fo,Fo',F5,F5) : retarded 
(Fo ,— Fq ,F5,— F5 ) : advanced 

- m 

(-Fo'fFJ,Fs,~FJ) : retarded 
(-F ',-F " ,F 5 ',F 5 " ) :advanced 

H<0, m>0 



FIG. 6: Fo and Fs as a function of m and fi. Fj and Fj 
are the real and imaginary part of a function Fj . At any pa- 
rameter point, we generally have at least two solutions, which 
correspond to the retarded Green function and advanced one. 




0.5 



f 






4" 










































0.05 ~~ 0.1 

m 




3. General fj, and m case 

For finite \i and m, both Fq and F5 are in general 
nonzero and we cannot solve eqs. (|38l4ip analytically. 
Accordingly, we have numerically solved the coupled 
equations with respect to a and b, so that Fq and F 5 
are derived in terms of \x and m. 

Before describing the numerical solution, let us first 
argue about the generic features of such solutions. No- 
tice first that ReG is an even function of both a and b, 
while ImG is an odd function of both a and b. Thus the 
following two should be degenerate at any given \i and 
to as the solutions of eqs. (|38ti4ip : 

(F^F^F^Fi'), (F ',-^',F 5 ',-F 5 "), (65) 

where Fj and Fj' are the real and imaginary part of Fj . 
Namely, these two solutions correspond to the retarded 
Green function and advanced one respectively. 

The above two solutions at given to and fi can be fur- 
ther extended into the other 3 quadrants, i.e. (— m, 
(to, — /i) and (— m, — 

(F , ±Fo , F 5 ', ±F 5 ")| m Al 

= (Fo , ±Fo , — F5 , T-F5 

= (— Fo 5 ±Fo ,F$ , T-F5 )\ m ,-n 
= (~Fq , ±Fq , — F5 , iFs )|_ m> _ M) 

where the upper sign corresponds to the retarded func- 
tion for any of these four regions by construction (Fig. [5]) . 
Observing this, please notice that both F§ and F§ vanish 
at m — 0, which is indeed the case with Sec.III-Al. Sim- 
ilarly, one can also see that Fq and F§ should vanish at 
fi = for any to. Both eqs. (|60H61I) and eq. actually 
observe this. 

These considerations are also consistent with the nu- 
merical solution. In Fig. [7] and [51 we demonstrated nu- 
merically how Fo and F5 behave as a function of \x and 
to (only for the first quadrant, m > and fi > 0), at 



FIG. 7: a — 0.48 ; (a): The contour plot of Fq as a function 
of n > and m > 0. The value of Fq decreases toward the 
dark region, and becomes zero at fi = 0. The contour interval 
is 1.2 x 10~ 3 . (b): The contour plot of Fq. The value of Fq 
decreases toward the darker region, and becomes zero at the 
yellow region. The contour interval is 0.6 x 10~ 3 . (c): The 
contour plot of F$. F^ decreases toward the darker region, and 
becomes zero at m — 0. The contour interval is 4.0 x 10 -3 . 
(d): The contour plot of Fg. Fl' increases toward the darker 
region and becomes zero at fi = 0, m = and the yellow 
region. The contour interval is —1.2 x 10 . 



specific values of a. Fig. [7] is for a < a c , while Fig. [5] is 
for a > a c . In the limit of — * 0, Fo and F5 in these 
two figures continuously connects with the two analytic 
solutions previously derived, i.e. eq. (|60ti61|) and eq. ([64]) 
respectively. 

We have also checked that, whenever Fq" = F§ = 0, F^ 
is always greater than Fq, i.e. a 2 — b 2 < 0. As such, the 
spectral function is identically zero, provided that both 
Fq and F§ vanish. Such a phase should be regarded as 
an incompressible phase having no bound states. On the 
one hand, when either Fq" or Fg' is finite, the spectral 
weight is finite and such a phase is compressible. 

By seeing whether Fq and F5' totally vanishes or not, 
we have determined the phase diagram in the fi-m-ot 
space. The phase boundaries between the compressible 
phase and the incompressible phase thus obtained are 
schematically drawn in Fig [9l while accurately specified 
in Fig. 1101 For a > a Cl we have a finite critical mass 
value, i.e. to c , below which a compressible phase extends 
even at n = (Fig. [1(c) and Fig. EDfe.f)). This critical 
value goes to zero, when a goes to ot c from above (Fig. 
UOf dlV For a < a c , we have a compressible region not 
in the fi = region anymore, but still in the nonzero 
fi region (Fig. H^b) and Fig. fTUJa-c)). The slope of the 
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FIG. 8: a = 0.52 ; (a): The contour plot of Fg as a function 
of fi > and m > 0. The value of Fq decreases toward the 
dark region, and becomes zero at /i = 0. The contour interval 
is 1.8 x 10~ 3 . (b): The contour plot of Fg'. The value of F," 
decreases toward the darker region, and becomes zero at the 
yellow region. The contour interval is 1.8 x 1CP 3 . (c): The 
contour plot of Fg. F5 decreases toward the darker region, and 
becomes zero at m = 0. The contour interval is 4.0 x 10~ 3 . 
(d): The contour plot of F$ . F$ increases toward the darker 
region and becomes zero at fi = 0, m = and the yellow 
region. The contour interval is —3.0 x 10 . 



phase boundary in a < a c given as follows; 

d[x c _ 1 — 2a 



dm. 



c \i_i c —m c —Q 



l + 2a' 



increases when the disorder strength decreases fFig HOT a- 
c)). 



IV. DIFFUSON AND QUANTUM 
CONDUCTIVITY CORRECTION 



In the previous section, we have derived the 1-point 
Green function within the self-consistent Born approx- 
imation. In the 3-d parameter space spanned by [i, m 
and a, we have observed that the topological insulator 
and an ordinary insulator are always intervened by the 
compressible phase (see the blue region in Fig. [5]) . The 
topological insulator supports a single 2 + 1 surface mass- 
less Dirac fermion on each boundary, while an ordinary 
insulator does not. As such, we expect that this inter- 
vening phase is composed by those wavefunctions which 
extend over an entire bulk (see section I for its reason). 

As the first step to understand the nature of this com- 
pressible phase especially for a < a c , we will calculate 
the series sum of the ladder- type diagram T d (q, to) (see 
Fig. [TTT a)) , only to see that the diffuson thus obtained 



FIG. 9: (a) A schematic phase diagram in the /i-m-a space. 
Either Fq or F§ always remains finite in the compressible 
phase (blue), while both of them become zero at the remain- 
ing parameter region (incompressible phase), which is further 
divided into an ordinary insulator (red) and the topological 
insulator (yellow), (b) A schematic phase diagram in [i-m 
plane for a < a c , and (c) that for a > a c . These correspond 
to the numerical results shown in Fig. \W\ 



consists of two quasi-degenerate low-energy modes; 



f d (q,tu)<x 



1 



u> + iDq 



2" 1 



w + iDq 2 + ir to ^ c 



(66) 



with r t ~p oc m? . The information of the charge diffusion 
is solely encoded into the first term, which thus always 
has the diffusion pole structure. On the other hand, the 
second term becomes massless only at m — (but generic 
/j,) , while it suffers from the finite infrared cutoff r t ~p for 
the finite m case. We will explicitly see that the second 
term is actually associated with the parity degree of free- 
dom, which, at m — 0, becomes a conserved quantity of 
our effective continuum model, i.e. eq. ([5]). 

When the hole line of T d (q,uj) is time-reversed, these 
two-mode features are translated into the backward scat- 
tering channel associated with the "fan" -type diagrams 
L> coop (fc + k',u) (see Fig. [TTJb)). Namely, for small uj 
and k + k', it is also dominated by two quasi-degenerate 
dominant contributions; 



U coop (k + k',u>) cx 



uj + iD(k + k') 2 
1 



Ul 



LO + tD(k + k') 2 +zr topo 



(67) 



One is obtained from the charge mode channel, i.e. Tf, 
with its hole-line time-reversed, while the other is from 
the parity mode channel, T d . In this section, we will fur- 
ther see that both of these two give rise to the same am- 
plitude of the anti-weak-localization (AWL) correction to 
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FIG. 10: Phase boundaries between the compressible phase 
and the incompressible (gapped) phase, in the /j, — m plane, 
at several values of a. (a) a = 0.0, (b) a — 0.1, (c) a = 0.4 
(d) a = a c = 0.5, (e) a = 0.6 (f) a = 0.7. 
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FIG. 11: (a) A series sum of the ladder-type diagrams 
V a -ysili w )i (b) A series sum of the "fan"-type diagrams 



the electric conductivity at m = 0. In the presence of the 
finite topological mass m, however, the second term in 
eq. (f6"7) becomes less dominant, because of the finite in- 
frared cut-off r t ~p . Namely, half of the AWL correction 
becomes ineffective on increasing m ("quantum correc- 
tion doubling" ) . 

Using these knowledges obtained in this section, we 
will propose in the next section the possible microscopic 
mechanism of how the bulk-critical region emerges be- 
tween the topological insulator and an ordinary insula- 
tor. 

This section is organized as follows. In Sec. IV-A, we 
will sum up the entire ladder type diagram in the particle- 
hole channel, using the 1-point Green function obtained 
in the self-consistent Born calculation; 



G^Qt; n,m) = F 1 - k v % + F 5 %, 
G A '- 1 (k;fi,m)=Ffl-k v %+F^ 5 , 



(68) 
(69) 



with Fq = ft+ir^ 1 and F$ oc m. Such a summand should 
contain those contributions which diverge at uj = and 
q = 0. We will identify this diverging contribution in 
the section IV-B, only to see that T d (q,uj) contain two 
gwosi-degenerate dominant contributions, as in eq. (|66[) . 
Explicit expressions for Tf , and r t ~p Q will be therefore 
given here. By calculating the parity-density correlation 
function, we will show in the section. IV-C that T d solely 
participates in the parity diffusion mode. Finally, the 
quantum conductivity corrections arising from these two 
terms are calculated in section. IV-D, based on the Kubo 
formula. 

Comparing eqs. (|53H54|) with eqs. (|57I59[) . notice also 
that the weak-localization (WL) calculation in this sec- 
tion becomes a controlled analysis only for the weak dis- 
order region, a < a c . Namely, the renormalized chemi- 
cal potential p, and the life-time r" 1 determined in sec- 
tion. Ill guarantee a sufficiently small l//2r around fi ~ 
only for this weak disorder region ; 



1 (a c — a) 2 1 / \ 



for at < ot c 



Q7T 



(70) 



For the strong disorder region, however, 1/fir readily di- 
verges around the zero energy region; 



2a 



for a > a c . (71) 



Thus, the NCA approximation employed in section. Ill 
and the corresponding WL calculation described below 
acquire the small coupling constant •=-, only for a < a c , 
but, for a > a c , they generally don't. Bearing in mind es- 
pecially this strong disorder region, we will demonstrate 
in the appendix. B the mode-mode coupling calculation, 
which is complementary to the weak-localization stud- 
ies described in this section. Without resorting to the 
Kubo formula, this theoretical framework gives us the 
gap equation for the diffusion constant, taking into ac- 
count the interference effects due to the Cooperon term. 
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The basic feature which we will observe in this section, 
such as the quantum correction doubling, are also con- 
firmed by this analysis (see Appendix. B). 



A. Ladder-type diagrams in the particle-hole 
channel 



For clarity of the following description, let us first de- 
fine a tensor composed by two 4x4 matrices; 

A- A.y x . 

The former 4x4 matrix A r is for the particle (retarded) 
line, while the other, i.e. A a , is for the hole (advance) 
line. Throughout this section, we distinguish this "x"- 
mark from the "<£>" -mark, latter of which separates the 
spin space and sublattice space. The product of two ten- 
sors is defined as follows; 



A ■ B = A r B r x B a A a . 



(72) 



Note here that the order of the product in the hole line is 
reversed, compared with that of the left hand side. Under 
this algebra, the series-sum of the ladder diagram in the 
particle- hole channel, i.e. Fig.fTlTab is just the inverse of 
the following tensor, 



{i-afl( q ,cj)} a 



5,7/3 



{i xi- tt VG s (0 + ) x G A (0-)} 



k 



a<5,7/3 



= S aS S 10 -a}^G^(O + )G^(O^ (73) 

k 

with (0±) = (fc±|,/i±^). Namely, the following identity 
can be readily checked; 

{i - «n(g,w)} aQi fii/3 r^ iJ>7ft ( 9 ,w) ee 8 aS 5 10 . (74) 

Notice also that we have already normalized the momen- 
tum by the ultraviolet cutoff A as in eqs. (|35H37p . so that 
a ee 27rAooA and ^ fc ee 2 J k 2 dk. This notation also 
holds true for eqs. (I78l92l93ll . 

For simplicity of the explanation, we will calculate the 
inverse of eq. (|73[) . with q taken to be zero from the be- 
ginning (see also the description around eq. (|99[) ). Such 
an inverse diverges if lo = 0, while it does not for gen- 
eral to. When q taken to be zero, the polarization part 
becomes simple; 

{i-an( g = o >w )} atfj7/J 

ee a l + aiTi + a 2 T 2 + a 3 T 3 + a 4 T 4 , (75) 



T\ = ^2 7m x %> f 2 = ix 75, 
T 3 = 75 x 1, T 4 ee 75 x 75. 



(76) 
(77) 



Namely, it is just a linear combination of the five tensors 
with their coefficients defined as follows; 



{ao, 01,02,03,04} 



'£ (h 



A, ' f , 

(k 2 - ^+ + Fi + W {FtJ 2 + (F*_)i) ' 

F 0+ F 5- 

(fc 2 - F 2 + + Fi + W ~ (^o-) 2 + (^5-) 2 ) ' 
F 5 +F£_ 



- ( fc2 _ f 2 + + Fi + w (F^r + (F 5 *_r) ■ 



^5+^5*- 



(fc 2 - F 2 + + F 2 + )(fc 2 - (F *_) 2 + (F*_) 2 ) 



.(78) 



The ± subscripts on Fi above mean that the argument 
of Fi((i) is shifted by ±|; 

F i ± = F i (jj l ±~). 

As such, Oq, Oi and 04 become real- valued at ui = 0, while 
02 and 03 become complex conjugate with each other; 



% 0^A\u,=0 - a 0,l,4|^ =0 , 02^=0 - O3| w =0- 



(79) 



One can evaluate the signs of the former three real- valued 
quantities, by noting that F^ is much smaller than the 
ultraviolet cutoff. They read, 



a o\u=o > 0, oi| w= o < 0, a4| w=0 < 0. 



(80) 



Notice also that all the integrands for 02, 03 and 04 con- 
tain F§, which is proportional to the topological mass 
(see sec. III). Thus these three quantities vanish when 
the topological mass m is zero; 



04 tx C(m 2 ), 02.3 cx 0(m). 



(81) 



The inverse of eq. (|75|) also becomes a linear combina- 
tion of a finite number of tensors composed of 7 matrices, 
because of the cyclic nature; 7172737475 = —70; 

f d (g = 0,c) 

= /3 1 + /?ifi + [3 2 T 2 + fof 3 + f3 4 f 4 
+f3 5 Sx + f3 6 S 2 + fofi ■ fx + fcf 2 ■ Sx 
+faf 3 ■ Sx + Pxafx ■ Sx + Pxx?4 ■ S 2 , (82) 

where two additional tensors are introduced in the fol- 
lowing way; 

Sx = 7i 72 x 7271 + 7273 x 7372 + 7371 x 737i, 
S 2 = 717273x737271- (83) 
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After lengthy algebra, one can express its 12 coefficients 
(3j in terms of those of eq. ([75]) as follows; 



A> 




— 3(5ao4 


-(5a 04 


— 3ao4 


— a 4 


Pi 




ai 


3ai 


ai 


3ai 


Pi 




35a 23 


&23 


3a 23 


023 


Pa 




-3<5a 23 


-<5a 23 


3a 23 


023 


Pi 


_ 1 


3<5ao4 


5a 4 


— 3oo4 


— ao4 


P-, 


= 8 




-<5a 4 


ao4 




06 




— 3ai 


3ai 


— 3ai 


3ai 


07 




-Ol 


— 3ai 


ai 


3ai 


Ps 




-^23 


&23 


~»23 


C123 


p9 




<5a 23 


-(5a 23 


-«23 


023 


010 




-<5a 4 


(5a 4 


ao4 


— a 4 




0u = - 


-3 ai 3 (/ 4 - 


1 + /r 1 )/i/ 2 ( 


-/s4 



/l 

h 
h 



(84) 



where ao4, <Sao4, a 23 , oa 23 and /i,2,3,4 are defined in terms 
of ao, ai, • • • and a^; 



<5a 23 


= a 2 


— a 3 , <5ao4 — ao — 04, 


(86) 


a 23 = 


= a 2 4 


- a 3 , ao4 = ao + 04, 


(87) 


fl = 




1 


(88) 


a? + 


(oa 4 + <5a 23 )(-<Sa 04 + <5a 23 ) ' 


h = 




1 


(89) 


a? + 


(—004 + ^23) (ao4 + 023) ' 


h = 




1 


(90) 


9a?- 


f (c5a 04 + <5a 23 )(-<5a 04 + <5a 23 ) ' 


h = 




1 


(91) 


9a?- 


f (-a 4 + »23)(ao4 + 023) ' 



B. Identification of the diffusion pole 

Using eqs. (|68H69I) . we have summed up the ladder- 
type diagram in the particle-hole channel, only to obtain 
eq. ([82]). The coefficients Pj appearing in eq. ([82]) are 
expressed in terms of Fq and F5, by way of eqs. (|84t - 
|9"T]) and eq. {75]). When the self-consistent Born (scB) 
solution is used for Fq and F5 , at least one of these fij is 
expected to have a diffusion pole structure. On the one 
hand, none of aj defined in eq. (|78p does not diverge at 
lu = 0. As such, some of fj~ should be zero at u> = 0. 
In this subsubsection, we will identify which fa diverges 
at small ui. This also determines the asymptotic tensor- 
form, of the diffuson in the small u) limit. 

To do this, let us first start from the self-consistent 
Born equations of F and F 5 , i.e. eqs. (|38H41[) . Or equiv- 
alently, begin with the following two; 



it? t? \ \ " (F0+F5) 



(F 2 -Fi) 



fj,-m + i5, (92) 
(j, + m + iS. (93) 



Then, subtracting eqs. (|92l93p by their complex conju- 
gates respectively, we can readily obtain 



O04 — «23 

3ai 



3ai 
ao4 + a 23 







-F><~ 







- |w=0 


FX 


+ F»_ 








,(94) 



where eq. (157]) and eq. ([75]) were used. This equation 
indicates that the determinant of the 2x2 matrix in the 
left hand side should be zero, provided that either Fq' or 
Fg' is non-zero. Any compressible phase having a finite 
density of state supports Fq' 2 + F5' 2 7^ 0. As such, any 
scB solution in the compressible phase always guarantees 
the following identity; 



{9a 2 - a; 



04 



*23J 



}|w=C 



0. 



(95) 



Since a\ defined in eq. ([78]) is negative definite at uj = 0, 
more accurately, eq. (]9"5]) should be replaced by; 



Observing eq. ([HI 
to the following, 



3ai| w=0 = -y/afa - a 32 |w=0 - ( 96 ) 
notice that this is actually identical 



0. 



(97) 



Namely, fa carries the diffusion pole. 

fa , / 2 and / 3 generally cannot have any pole structure 
for the small u> region. To see this explicitly, note first 
that, when generalized into the finite q case, fa takes the 
following asymptotic form; 



AM 



1 1 

aoT iuj 

fa(q,u 



1 



1 



oqt iuj — Do 



(98) 



In the right hand side, we have replaced ito by iui — D$q 2 
with the bare diffusion constant Dq. By retaining the 
subleading contribution in small q appearing in eq. (]73p . 
one can explicitly calculate its leading order expression 
in the large fir limit; 



Do 



1 a c 



6 a c + a 



(99) 



which is positive definite for a < a c . Similarly, we can 
obtain the asymptotic form of fa, fa and fa at to, q ~ 0; 



1 



1 



fa 

h 
h 

T-f 1 and t^ 1 above are positive definite; 



aor iu — 


D ig 2 


-rr 1 


1 


1 




aor ilu — 


D 2 q 2 




1 


1 




aor iuj — 


Dq 2 - 


-1 

^topo 



(100) 
(101) 
(102) 



'top 



o>0, 



8a? 
a 



>0, 



(103) 
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while r t0 p being positive semi-definite; 

4(— a a 4 + a 2 a 3 



topo 



(J 



> 0. (104) 



|o>=0 



The two inequalities in eq. (| 1 03|) and eq. (] 1 04[) are in- 
deed supported by eq. ([75]) and eqs. I|79ll80p respectively. 
These expressions indicate that fa, fa and fa always ex- 
perience the infrared cutoff for the low-energy and long 
wavelength region. 

Comparing eq. (|104p with eq. (fgTj) . notice also that 
r t Qp reduces to zero in the massless case, T t ~p cx m 2 , 
since a 2 , 03 and a 4 being zero. As such, fa acquires a 
same diffusion pole as fa does in the absence of the topo- 
logical mass. Meanwhile fa and fa always suffer from 
the (relatively large) finite infrared cutoff t^ 1 , irrespec- 
tively of the topological mass term. Thus, we will retain 
in eq. (|84|) only those terms proportional to fa and / 4 . 
Based on the same spirit, we will also replace eq. (|8"5j) by 
its leading order contribution in small u> and q; 



/ 9 u ~3ai(-/ 3 + / 4 ), 



(105) 



and use the following asymptotic expressions for 
{flo, ■■■ , An}', 



{0o, 04, 05, /?io} 

^0041^=0 



/ 3 {-l, 1,-1,1} 



ao4| w= o 



Mi, 1,1,1}, 



{0i,0 6 ,0 7 ,0n} 



^Mi,i,-i,-i} + ^Mi ) i,M}, 



{02, 



/ 3 {1, -1,1,-1} 



8 

a 23|^=Q 

8 



/ 4 {i, 1,1,1}. 



With these equations, the asymptotic form of the dif- 
fuson in small ui and q will be derived out of eq. (|82[) . It 
consists of the two gwasi-degenerate dominant contribu- 
tions; 



(106) 



where the two a;, g- free tensors are given as follows; 



1 1 



-a,Q4 3ai 023 



j=0 



1 2 



— <5a 04 3ai <5a 23 



|w=0 



(1 + T 4 ) • (1 + 5i) 

(i + f 4 )-(f 1 + ^ 2 ) 
(f 2 + f 3 )-(i + 5 1 ) 

(i - f 4 ) ■ (i + &) 
(i - f 4 ) • (fi + 5 2 ) 
(f 2 - f 3 ) • (i + 5i) 



(107) 



(108) 



C. parity diffusion mode 

To capture the physical meanings of the two members 
in eq. (|106p . notice first that, in the absence of the topo- 
logical mass m, our hamiltonian, i.e. eq. ([5]) with the 
chemical-potential type disorder, becomes invariant un- 
der the following U(l) transformation; 



71. (109) 



irrespectively of the strength of the disorder. As a re- 
sult, each ensemble at m — acquires another conserved 
charge, ip' W 745 "0(0, which is the parity density degree 
of freedom (see Table. II). Observing this U(l) symmetry, 
we can then expect that the diffuson T d (q,ui) calculated 
above should consist of two diffusive modes at m = 0: 
One describes the usual charge diffusion and the other is 
for the diffusion of this parity density. These two physi- 
cal modes actually correspond to the first term and the 
second term in eq. (|106p respectively. In fact, the par- 
ity density becomes non-conserved in the presence of fi- 
nite m, which is consistent with the finite infrared cut-off 
r t Qp cx m 2 appearing only in fa (see eqs. (11021981) ). 

To uphold this consideration more directly, one can 
also calculate the density correlation function and parity 
density correlation function at m = 0; 

(f>o(q,u) = ^2 ®a0,P<x(k,k';q ) Uj), 
k,k',a,/3 



^45(<?,w) = ^2[% 5 ] pa ^ a s,-yi3(k,k';q,uj)['y4 5 ] 



67' 



where <3E>(fc, k'\ q,ui) stands for the response function (see 
eq. (|B2[) for its definition). By noting that this response 
function for small q and ui is dominated by the diffuson; 

$ a s,y0(k,k';q,u) ~ -— G* ai {k+,n + )Gj ll3 {k-,ii-) 

x[r d (9,^)] ai , 1)7l/3l Gf l5 «,M+)^ 7l (^, M -), (HO) 

one can explicitly see that the two terms appearing in 
eq. (|106p actually contribute the density correlation and 
parity density correlation separately; 

(Mq,"),<f>' i5 (q,oj)) ^-^x (/ 4 ,/ 3 ). (111) 



D. 



Cooperon and the quantum conductivity 
correction 



When the hole lines being time-reversed, eq. (|106p will 
be transcribed into the two quasi-degenerate dominant 
contributions to the series sum of the "fan" -type dia- 
grams (see Fig. [TTJb)); 



C/ c °°P(fc + k', w) = ^-{faUt + faU c 2 ) 



\q— >fe+fc' ' 



(112) 
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where the u>, q-free tensors 
eqs. (I107fll08|) respectively; 



C/f 2 are derived out of 



Ul = 



-ao4 3ai a 23 



\uj=0 



(i + f 4 )-(i- 

(i - f 4 ) • (-f a 

(f 2 + f 3 ) • (i - 



Si) 
f S 2 ) 

Si) 



,(113) 



— <5a 4 3ai <5a 23 



j=0 



(1 - T 4 ) • (1 - Si) 
(l + f 4 )-(-fi+S 2 ) (|114) 

(f 2 - f 3 ) • (i - Si) 



Substituting these two Cooperon terms into the current- 
current correlation function, we can explicitly show that 
the two members in eq. (I112|) lead the same magnitude of 
the anti-weak-localization (AWL) behaviour at the criti- 
cal point (m = 0); 



E -f 1 E^ A ( fc ' /*) • 7i • M)}^r, Q/3 , 7 ,{G fl (-fc, ^ • 7 i • c A (-k, ^ 

L- 1 <\k+k'\<l- 1 k 

= E If E^( fc > • 7i ■ m)}^!^^^ M) • 7i • /x)}^ 

= f E (fc + fc/p E^ A ( fc ' ^) • 7i • G fl (fc, M )} 5a {i - Si - fi + S 2 } a ^ s {G R (-k, M ) • 7i • G A (-fc, /x)}^ 

^cxfr'-r 1 ) (115) 

i 



with c = 167r • (1 + — ) being positive definite, (we used 
0-2 = «3 = »4 = and — 3ai = ao in eqs. (|113H114[) ). 
When the finite topological mass to is introduced, how- 
ever, the second member of eq. (|112[) becomes suppressed, 
since the infrared divergence of becomes truncated by 
finite r t ~p . As a result, one half of the AWL correction 
becomes ineffective in the presence of finite to ( "quantum 
correction doubling" ) . 

V. DISCUSSION 
A. Summary of our findings 

In this paper, we have studied the effects of the time- 
reversal invariant disorder on the quantum spin Hall sys- 
te m 26 i 27 i 28 i 29 i 30 . We have especially focused on the quan- 
tum critical point (QCP) which intervenes the 3-d topo- 
logical insulator (TI) and an 3-d ordinary insulator. The 
topological insulator supports a single 2 + 1 massless sur- 
face Dirac fermion for each boundary, while an ordinary 
insulator does not have any. As such, the bulk wave- 
function in those parameter regions (or point) which in- 
tervene these two insulating phase should be extended, 
so as to mediate two opposite surfaces. Such extended 
bulk states are stable against T-invariant disorders, as 
far as each surface state in the TI phase is stable. In 
fact, Nomura et al^ and Bardarson et al£^ have recently 
calculated the (3 function numerically, and demonstrated 
that the single-copy of the 2 + 1 massless Dirac fermion 
is topologically stable against the T-invariant disorders. 
This observation strongly indicates that that there always 
exists delocalized (bulk-critical) region between the 3-d 



topological insulator phase and an 3-d ordinary insulator 
phase. 

To uncover the nature of this peculiar quantum criti- 
cal point (or region), we have studied the disorder effect 
on its minimal model, i.e. the 1-copy of the 3 + 1 Dirac 
fermion. As a basis for this, we first studied in the sec- 
tion III how the chemical potential type disorder brings 
about a finite life time of the zero-energy wavefunction. 
We then observed that there exists a certain critical dis- 
order strength above which the DOS at the zero-energy 
becomes finite (see eq. ([BT])). 

When the finite topological mass is introduced, a sys- 
tem eventually enters either the TI or an ordinary insu- 
lator, depending on the sign of the topological mass. In 
Sec. Ill, we studied how this topological mass are renor- 
malized by the chemical-potential-type disorder within 
the self-consistent Born approximation. By doing this, 
we have determined the phase boundary between the 
compressible phase and the gapped phase (see Figs. 19110)) . 

To further infer the low-energy structure in this com- 
pressible phase, we have derived in the section IV the 
diffuson, Cooperon and the weak localization (WL) cor- 
rection to the electric conductivity. We then observed 
that the charge diffusion mode and parity diffusion mode 
dominant the diffuson (see eq. (|106|) '); the charge chan- 
nel always carries the diffusion pole structure, while the 
parity density channel becomes massless only in the ab- 
sence of the topological mass. In the presence of the finite 
topological mass to, it generally suffers from the infrared 
cut-off t~* cx to 2 . 

Corresponding to this feature in the diffuson, the 
Cooperon is also composed of two quasi-degenerate domi- 
nant contributions (see eqs. (|112j|114j0 . In the zero topo- 
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logical mass limit, these two contributions bring about 
the same magnitude of the anti- weak-localization (AWL) 
correction with each other. When the finite topologi- 
cal mass m is introduced, however, that from the parity 
density channel becomes truncated by the finite infrared 
cutoff r t ~p . As such, one half of the AWL correction be- 
comes ineffective. As a result, on increasing m, the AWL 
correction exhibits a crossover into one half of its original 
value ( "quantum correction doubling" ) . 



B. "Levitation and pair annihilation" phenomena 

Let us discuss open issues in the 3-d Zi QSH system in 
the view point of our findings. As a tightly related topic 
to the stability of the QCP, the levitation and pair annihi- 
lations phenomena of the extended states^ were recently 
observed in the 2-d Zi quantum spin Hall systems by On- 
oda et al*& They numerically studied the disorder effect 
on the Kane-Mele model^ on the honeycomb lattice. In 
the clean case the system is set to be in the QSHI phase; 
namely, the spectrum consists of two bands, and there 
is a gap between them. When the system is disordered, 
some states far from the band centers become localized, 
while there are energy regions of delocalized states, lo- 
cated at the centers of the upper (empty) band and lower 
(filled) band. What Onoda et al. 26 have found is that 
each of these two does not disappear by itself, when the 
disorder strength is increased. Instead, when the disor- 
der becomes much stronger than the disorder strength 
for the localization in an ordinary insulator, these two 
merge into one bundle of extended states energetically, 
and annihilate in pair (see Fig. fT2Ta)h 

To argue this phenomena more generally, consider the 
3-dimensional parameter space spanned by the topolog- 
ical mass term m, chemical potential fi and disorder 
strength a. From the surface-state arguments described 
in the introduction, two insulating phases having differ- 
ent types of edge (surface) states, i.e. the topological 
insulator and an ordinary insulator, should be discon- 
nected by the delocalized (bulk-critical) region. Then, 
when a finite topological mass term m changes its sign 
from positive (topological insulator side) to negative (or- 
dinary insulator side), we should also expect that a sim- 
ilar levitation and pair annihilation phenomena occurs. 
Namely, when a system transits from the topological in- 
sulator side to the ordinary insulator side, the region of 
extended states in the upper band and that in the lower 
band always merge and annihilate with each other (see 
Fig. fT^T b)). Combining this picture with the Onoda's 
numerical observation^, one can then expect that the 
delocalized (bulk-critical) region constitute a surface in 
the 3-d parameter space spanned by fi, m and a, only to 
isolate the topological insulator phase from an ordinary 
insulator phase (see Fig. IT5]l . 
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FIG. 12: (a); Onoda's phase diagram— in the a-fj, plane with 
m > (topological insulator side) (b); A schematic phase 
diagram in the m-fi plane, which is expected from the sur- 
face (edge) state's argument described in the section. I. (c); 
A schematic phase diagram in the a-/i plane of the ordinary 
insulator side ( m < ). In (a-c), we have two delocal- 
ized energy regions (blue filled regions), which locate at the 
center of the upper band and the lower band. In (a), these 
two delocalized regions eventually merge and annihilate in 
pair, when a increases. As a result, the topological insulator 
and ordinary insulator are always disconnected by the bulk- 
critical (delocalized) region. In (c), however, two delocalized 
regions registered at the upper band and the lower band an- 
nihilate without merging each other, when a increases. Thus, 
all the insulating regions appearing in (c) are adiabatically 
connected from one point to others. In (b), two delocalized 
regions merge and annihilate with each other, when the topo- 
logical mass m changes its sign from positive to negative. As 
a result, the topological insulator and ordinary insulator are 
again disconnected from each other by the bulk-critical (de- 
localized) region, as in (a). 
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FIG. 13: A Schematic phase diagram in the a-fi-m space. The 
vertical axis is the disorder strength a, while the horizontal 
plane is spanned by the chemical potential fi and the topo- 
logical mass m. The delocalized region (blue filled region) 
constitutes a surface in this 3-dimensional parameter space, 
so that an ordinary insulator phase and the topological in- 
sulator phase are adiabatically disconnected from each other. 
Namely, one cannot move from one phase to the other, with- 
out crossing the delocalized region, i.e. bulk-critical region. 
The phase diagram for the constant positive m (topological 
insulator side) and that for the constant negative m (ordinary 
insulator side) are separately described in Fig. I12f a) and (c) 
respectively. The phase diagram for the constant a corre- 
sponds to Fig. fT^T b). 



C. Possible microscopic scenario 

Generally speaking, one has to go beyond our mean- 
field treatment of disorder in order to study the behav- 
iors of mobility edges. However, we can still speculate 
the microscopic picture of the "levitation and pair an- 
nihilation" phenomena discussed above, in terms of the 
"quantum correction doubling" found in this paper. 

We expect that the intervening bulk-critical region 
(blue filled region in Figs. [TWa-b) and Fig. \T3§ corre- 
sponds to the T t Qp Q = region. Namely, when a sys- 
tem transits from the topological insulator to an ordi- 
nary insulator, we surmise that one of the high-energy 
modes, i.e. parity diffusion mode appearing in eq. M06}) , 
becomes massless once, only to guarantee the existence 
of the bulk-critical region between these two insulating 
phases. This conjecture naturally leads to the following 
microscopic scenario of the "levitation and pair annihila- 
tion" phenomena. 

Suppose that the T-symmetric disorder is introduced 
in the topological insulator. We assume that such dis- 
order potential is strong enough to make the system lo- 
calized. But it is not strong enough to make the upper 
(empty) band and low (occupied) band mixed with each 
other Namely, a system locates in the topological insula- 
tor side of Fig.fTWa). so that it can be adiabatically con- 
nected into the topological insulator phase in the clean 



limit. In such localized phase, we expect that a parity dif- 
fusion mode always exists in the high energy region and 
is protected by the infrared cut-off T t ~p , while the charge 
diffusion mode disappear because of the relatively strong 
disorders (see Fig. [T4ja)). Starting from such localized 
phase, decrease the topological mass term (or further in- 
crease the disorder strength). Then, this infrared cutoff 
r t Qp associated with the parity diffusion mode decreases 
gradually, only to be renormalized to be zero at the tran- 
sition point (see Fig. [TWb)). Namely, at this transition 
point, the parity diffusion mode becomes massless. As 
a result, the Cooperon term corresponding to this par- 
ity diffusion mode, i.e. eq. (|114[) . becomes effective and 
brings about the positive quantum interference effect on 
the back-scattering processes, in a same way as in the 
section. IVD. Because of this positive quantum interfer- 
ence effect, which emerges only when the parity diffusion 
mode becomes massless, the charge diffusion constant re- 
covers at around r t ~p ~ 0, even in the presence of the 
relatively strong disorder (see the red line in Fig. [TWb)). 

However, when one further decreases the topological 
mass (or increases the disorder strength), the infrared 
cutoff r t ~p becomes finite again. As a result, this posi- 
tive quantum interference effect due to the massless par- 
ity diffusion mode becomes ineffective again, only to lead 
a system into an insulating phase (see Fig. [TWc)). This 
insulating phase is now adiabatically connected to an or- 
dinary insulator in the clean limit. 

To uphold this microscopic picture, we need to con- 
sider several ingredients missing in our approach 35 . We 
will enumerate them in the following. As indicated in 
Figs. [12]fT3l the pair-annihilation occurs only in the topo- 
logical insulator side. Namely, the phase diagram is 
asymmetric with respect to the sign change of the topo- 
logical mass term. On the other hand, all the findings in 
this paper are symmetric with respect to the sign change 
of this mass term. This is obviously because our starting 
model is the effective continuum model, describing only 
the local structure around a certain fc-point. On the other 
hand, the Zi topological number is determined from the 
global information of the Bloch wavefunctions' phase in 
the k spaced. Therefore, in such an effective continuum 
model one cannot determine whether the topological in- 
sulator by itself corresponds to the to > phase or the 
to < phase. Instead, it simply dictates that one of 
these two should be the topological insulator, and the 
other is an ordinary insulator. As such, to describe the 
asymmetric behavior of the mobility edge as in Fig. [T3J 
we clearly have to deal with lattice models. 

In the present work we treated disorder in the mean- 
field level, considering only the Cooperon correction. To 
verify the aforementioned scenario, we thus also need 
to deal with interactions among the various low-energy 
modes, beyond the mean-field treatment. In such sit- 
uations, the inter-mode interaction between the quasi- 
degenerate Goldstone modes found in the section IV cer- 
tainly plays an important role in the "levitation and pair 
annihilation" phenomena. 
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FIG. 14: (a); The low-energy spectrum in the topological in- 
sulator side contains two relevant diffusion modes. One is the 
usual charge diffusion mode, which disappears in the presence 
of relatively strong disorders. The other is the parity diffusion 
mode, which is protected by the infrared cut-off T t ~p from the 
disorders, (b); When a system transits from the topological 
insulator phase to the ordinary insulator phase, the parity 
diffusion mode becomes massless. Namely, when one further 
increases the disorder, starting from (a), r t ~ po becomes renor- 
malized by the disorder, only to reach zero. As a result, the 
Cooperon term corresponding to this parity diffusion mode 
becomes effective and induces the positive quantum interfer- 
ence effect on the backward scattering process of the charge- 
degrees of freedom. Because of this, the charge diffusion mode 
recovers at r t ~p Q = 0. (c); The low-energy spectrum in the 
ordinary insulator side. 



APPENDIX A: EFFECTS OF GENERIC 
TIME-REVERSAL INVARIANT DISORDERS 

In this paper, we have restricted ourselves to the chem- 
ical potential type disorder for simplicity. However, there 
exist in general several other types of T-invariant disor- 
der potentials, as described in section II (see eq. (0)). 
We basically expect that these additional time-reversal 
invariant disorders will not change our results drasti- 
cally. To uphold this expectation, we study in this ap- 
pendix how our self-consistent Born solution is modified 
in the presence of generic time-reversal invariant disor- 
ders, focusing on the zero-energy wavefunction at the 
critical point. As a result, we will obtain the follow- 
ing facts, which support this expectation. One is that, 
when only the diagonal correlations, Ajj, are present, 
our solutions derived in section III do not change at all 
(see eas.f IA16IA4IA8)) ). When the off-diagonal correla- 
tion such as A05 is introduced, F§ acquires a finite imag- 
inary part, i.e. F§ ^ 0, even at the zero-energy state 
at the critical point (see eqs. (|A16IA4IA8|) h However, 
provided that A05 is not so large in comparison with the 
diagonal correlation such as A o, A 55 and etc. the effect 
of the non-zero F§ is expected to be negligible. 

The generic T-reversal invariant disorders bring about 
the coupling between F5 and Fq more explicitly in the 
self-consistent Born (scB) equation. Namely, instead of 
eqs. (I38H39[) . our scB equation reads; 

f ** £ +F °- 5 BF ; 2 =fo-F , (Al) 
Ja<\k\<i Fti-Y,„=iFjl 

i d3k T^~ F L = h - F5 > (A2) 

Jo<\k\<i Jo-E^i*? 

where only the following three parameters are the rele- 
vant model parameters; 

A ± = {A o + A 55 ± A «} A < ( A3 ) 

je{15,-,45} 

B = 2A 05 A. (A4) 

The coefficients of 71,2,3,4 in the 1-point Green function, 
on the other hand, are again free from renormalization; 
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iT,2,3 = A,2,3 = -h,2,3, F 4 = f 4 = 0. (A5) 

In terms of G defined in eqs. (|40I41[) . We can rewrite 
eqs. (|AltiA2|) more transparently; 

2ir(A+F Q - BF 5 ) ■ G = f - F , (M) 
2K(BF -A-F b )-G = h-F b . 

When it comes to the zero-energy wavefunction at the 
critical point, i.e. /o = /s = 0, this coupled equation 
could be "diagonalized" ; 

(l-77 ff G)(A CT F -F 5 ) =0, (A7) 
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with a — ±. rj a and A CT are denned as follows 
1 



A ± = -(A s ±^Af 



B' 2 



v± = -(A Q ± V /A2- J B2) ) 
with positive definite A s and A a ; 



(A8) 
(A9) 



A, = \{A + 
Aa = \{A + 



A) = A 00 A + A 55 A 



j£{15,-,45} 



A W A. 



Observing eqs. (fT7|) . note also that A s defined above is 
always greater than \B\ defined in eq. (|A4[I ; 



A 2 ~B 2 > 0. (A10) 

Eq. (|A7|) with ct = ± can be trivially satisfied by Fq = 
fg = 0. In what follows, we will enumerate all possible 
non-trivial solutions of this coupled equation. Let us first 
write down the real part and imaginary part of Eq. (|A7j) 
for both a = ±, separately. Noting that A± and rj± are 
real- valued, we have the following for a = +, 



1 — 7/+ReG 77+ ImG 
— 77+ImG 1 — 7/ + ReG 

For a = — , we have 



l-?7_ReG ?7_ImG 
-))JmG 1 - ?7_ReG 



A_itf-F 5 " 



0. (All) 



0. (A12) 



Observing eq. (|A10[) . notice that A + ^ A_ in general. As 
such, (F ,F 5 ) cannot satisfy F 5 = \-F and F 5 = A + F 
simultaneously. Thus, when i*g — \-Fq is adopted, the 
determinant of the 2x2 matrix in eq. (|A1 1[) should be 
zero; 



1 — ?/ + ReG ?7 + ImG 
— ?y + ImG 1 — ?7+ReG 



= 0, 



(A13) 



or equivalently 



1 = 7^+RcG, ImG = 0. 

On the other hand, when F5 = X+Fo holds true, we have 
the following in turn, 

1 = ?y_RcG, ImG = 0. 

We thus have the only two possible non-trivial solutions; 

f (Bi) : F 5 = A_F , 1 = ry+RcG and ImG = 0, 
\ (Bii) : F5 = A+F , 1 = ?7_ReG and ImG = 0. 

In either cases, ImG = readily leads us to a = first. 
The reasoning of this was already described in section 
IIIA1. When a = 0, the real part of the function G 




0.0 
B/As 



B/As 



FIG. 15: The phase diagram of the scB solution in the pres- 
ence of generic time-reversal invariant disorders. The region- 
(a) includes the "compressible phase" argued in the section III 
and IV, i.e. a > a c and B — and A a . The region-(b) ap- 
pears only when A a = X/j={i5 25 45} ^ii -* Note also 
that A a , A 3 > and B < A s because of eq. (|17|l (see also 
eq. <[XT0])h 



becomes simplified; ReG| a=0 = — 2 + 26ArcTan[6~ 1 ] (see 
eq. ([i5|) h Thus, above two solutions will be transcribed 
into the following two; 



6ArcTan[6 _1 ] = 



1 + 2t?, 
2ry T 



(A14) 



with F5 = X±Fq respectively. 

Since the left hand side of eq. (|A14[) is positive semi- 
definite, we have the following two parameter region sup- 
porting non-trivial solutions; 

J (a) : V+ < -5 < 
\ (b) : TJ+ < //_ < -\. 

We also used ?/_ > which is trivially supported by 
eq. (|A10|) . These two parameter regions are depicted in 
Fig- EH where the region- (a) actually includes the "com- 
pressible phase" argued in the sections III and IV, i.e. 
a > a c and B = A a = 0. 

In this region-(a), only the type-(Bi) solution becomes 
possible; 



F 5 = A_F , (a,6)~ (0,. 



27T774 



(A15) 



Under F^—F 2 = (a+ib) 2 , this is identical to the following 
solution; 



(F ,F 5 )=±i 



\b\ 



1-X 2 



=(1,A_). 



(A16) 



This solution comprises continuously with the physical 
scB solution described in the section III. Namely, when 
A a taken to be zero, eq. (|A16j) precisely reduces to 
eq. (EU). 
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When it comes to the region-(fo), type-(Bii) also be- 
comes a possible solution; 

F 5 = \ + F , (a, 6)^(0 /+^- ), (A17) 

namely, 

{F ,F B ) = ^-=B=(\+ 1 ,1). (A18) 

vV - 1 

In the absence of finite B and A a , however, this solution 
is continued into eq. (|49|) . Thus, this can never hold true 
therein. Because of this, we judge the type-(Bii) solution 
to be unphysical. 



has the small coupling constant 1 / (fir) only for the 
weakly disordered region, i.e. a < a c , while it be- 
comes an uncontrolled analysis for a > a c . Bearing in 
mind this strong disorder region, we will employ in this 
appendix more phenomcnological calculations, based on 
the mode-mode coupling theory 23 . Without resorting to 
the Kubo formula, this theoretical framework gives us a 
mean-field equation for the diffusion constant D, where 
the quantum correction due to the Cooperon term are 
taken into account as in the standard weak-localization 
calculation^ 3 -. The final results of this appendix. B such 
as eqs. (|B54IB61[) indicate that this quantum correction 
becomes doubled, when the topological mass term m are 
fined tuned to be zero. 



APPENDIX B: MODE-MODE COUPLING 
THEORY 

The weak-localization calculation (and also the self- 
consistent Born calculation) described in the section. IV 



The calculation consists of two steps. The first step 
begins with the Bethe-Salpeter (BS) equation for the re- 
sponse function & a 5,-f/3(k, k'\ q, lu); 



G« Ql (fe+,M+)G^(fc-,M-){ - ^-<W<Wfc,fc' +^Ul P J* ni(3l (k,k 1 ;q,uj)<5>s 1 S n - <1 (ki,k'-,q,u J )y (Bl) 
*ai, 7 /j(*, k'- q, u) ee -^-(G* 5 (k+, k' + , fi+)G^(k'_, fc-,M-)> imp , (B2) 



with k± ee k ± I and fi± = // ± -| . Out of this equa- 
tion, we first derive the linearized equations of motion 
(EOM's) for the density relaxation function 4>o(q, w), cur- 
rent relaxation function tfrj(q,u>) and relaxation functions 
associated with other internal degrees of freedom; 

<M<?,W) EE ^ [7o] /3Q $a 7 ,7/3(M';<? I w) I ( B 3) 

4>MiU) ee Y %[%] /3a ^a J ^0(k,k';q,Uj)(BA) 

k 7 k' a, 0,7 

with q normalized to be a unit vector. Since the EOM's 
thus derived are linearized, one can solve them for these 
relaxation functions, only to obtain their asymptotic ex- 
pressions for the small q, o>; 



<t>a{q,u) 



1 



U) + iDq 



U) + iDq 2 ' 



(B5) 



To be more specific, the diffusion constant D appearing 
in the denominators will be expressed in terms of the re- 
laxation kernels M a .b(q, latter of which are defined by 
the two-particle irreducible (2PIR) vertex function [/ 2PIR 



(step-(i)); 



D EE i 



M 



5j,5j 



M ,.,Mr,, r,, - MsjjMjfii J {q ^ 



(B6) 



M a , b (q,ui) =2ir 1 S ab + — —x 

Z 7T 



Y ba(k; q, co)} 0a U^(k, k'; q, u) [ 7fc R (fc'; q, W )] ^(B7) 
k,k' 

(see also eqs. (|B3HIB32I) for the definitions of 7a 3 ). 

The 2PIR vertex function U 2P1R (k,k';q,uj) in disor- 
dered media is usually dominated by the Cooperon at 
small uj and k + k' . The Cooperon is the series-sum of 
the ladder-type diagrams in the particle-particle channel, 
which is therefore obtained from the diffuson with the 
hole-line time-reversed. The diffuson is in turn responsi- 
ble for the diffusion pole in the relaxation functions, i.e. 
the denominators in eq. (|B5[) . As such, in the presence of 
the T-symmetry, the asymptotic form of the Cooperon 
at small to and k + k' should be characterized by the 
same diffusion constant as that in eq. (|B5j) . Based on this 
spirit, we will replace the 2PIR vertex function in eq. (1B7|) 
by this asymptotic form of the Cooperon. Through this 
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approximation, eqs. (|B6tiB7|) constitute a self-consistent 
equation for the diffusion constant D (step-(ii)). 

As we have seen explicitly in section. IV, the diffuson 
consists of the charge diffusion mode and parity diffusion 
mode; 



T d (q,Lo)cx 



1 



iDq- 



1 



iDq 2 + ir t 



— r2,(B8) 



topo 



with the positive semi-definite T t ~p proportional to to 2 
(see eqs. (|104I107H108[) V Namely, the second term, i.e. 
the parity diffusion mode, generally suffers from the finite 
infrared cut-off in the presence of the topological mass, 
while both of these two equally dominates the low-energy 
region for to ~ 0; 



,(ff- 



t d (q,u)<x{ ZJTm ^_ 

{ ZJ+lDq 1 1 



for DV 
for DV 



> T, 



topo ) 

.-1 

topo 5 



(B9) 



(see eqs. (|102I98[) and eq. (|106p ). In the presence of the 
T-symmetry, this crossover behaviour will be transcribed 
onto the Cooperon term; the backward scattering process 
originated from the parity diffusion mode becomes inef- 
fective, in the presence of the relatively large topological 
mass to; 



U co °P(k + k',cj) oc 



ui+iD(k+k')'^ 



for Dl- 2 <T t - l , 
for £»r 2 >r t ^ , 



(see eq. ([XI2])). 

Corresponding to these two-mode features, we will de- 
rive in this appendix the two limiting gap equations; 
one is valid for Dl~ 2 > r t ~p , while the other is for 
Dl~ 2 < T t ~p Q . This appendix is organized as follows. 
The appendix. BI is devoted for the step-(i), in which 
the linearized coupled EOM's for the relaxation functions 
and eqs. (|B6flB7|) will be derived. Using eqs (|B64IB7j> . 
we will derive in the appendix. B2 the gap equations 
for the two- limiting cases (step-(ii)). By solving these 
gap equations, we will finally see how the diffusion con- 
stant for a > a c behaves as a function of ft and to (see 
eqs. (|B54IB6ip and Fig. EH]). 



1. Coupled EOM's for relaxation functions 

The EOM's derived henceforth are linearized with 
respect to the relaxation functions (unknown quanti- 
ties). Namely, the mode-mode interactions among vari- 
ous bosonic degrees of freedom (density, current and so 
on) will be represented by the 'mean-field' induced by the 
corresponding relaxation functions. This mean-field for 
the relaxation function is analogous to the self-energy for 
a I-point Green function, so that it is often dubbed as the 
relaxation kernel 23 . Being linearized, such EOM's can be 



easily solved, only to let us express relaxation functions 
in terms of the relaxation kernel. 

These linearized EOM's also have to be closed with re- 
spect to a set of unknown relaxation functions. Consider, 
for example, the EOM of the density relaxation function, 
which is nothing but the continuity equation. This equa- 
tion contains the current relaxation function. Accord- 
ingly, to make coupled EOM's to be closed, we further 
need the EOM for this current relaxation function, i.e. 
constitutive equation. The constitutive equation usually 
involves interactions between the current and other de- 
grees of freedom (DOF's) such as the spin density, sub- 
lattice density and so forth. As such, we further need to 
derive the EOM's of the relaxation functions associated 
with these internal DOF's. In this way, we need to make 
our entire coupled EOM's to be closed with respect to a 
set of unknown relaxation functions. 

Let us begin with the continuity equation. Apply the 
following differential operator from the left hand side of 
the Bethc-Salpeter equation eq. (|B 1[) : 

5 G~ 1 (k;q,uj) = G R '- 1 (k+, M +) - G A <-\k-, //_) 
= wl - 9 a7a - t R {pi+) + S a (m-)- 

Taking the summation over repeated band indices, we 
then have; 

[ujI - qpTYp, - Z R + S A ] /3a $ a7 , 7 ^(A:, k'; q, u) 

= -[G R (k + ,n + ) - G A (k.,^)] piai { - ^-5 aiPl 5 Kk , 

+ ^Xtia fc i; ?> ^$5x7,77! w) jaio) 

Under the integrals over k and fc', the vertex function 
and the self-energy in eq. (|B10|) set off each other; 

uj(/)o(q,uj) - q<pj(q,u) 

hi £ i d *-yK^+) - G A 7 (k'_,^)}. (Bll) 



k' ,7 



Namely, we used the following Ward identity; 

[s fl (fc + , M+ )-s A (fc_, / ^)] /3a 

= £ *<W ( fc '; 1, u)U 2 Jl%, (k\ k; q, uj), 



with 6G(k;q,uj) = G R (k+,^ + ) - G A (fc_,^_). 

Recall that we are interested in the relaxation func- 
tions for sufficiently low-energy and long wave-length re- 
gion; only to derive their diffusion pole structure. Thus, 
regarding u> and q as sufficiently small quantities, we can 
replace the right hand side of eq. (jBll [) by the spectral 
function; 

w0o(g,w)-g0j-(g,w) = i4o + O(g,w), (B12) 
where \Aq\ stands for the density of state at e = fi; 

A = -^rir[5G(fc;0,0)] ~ -16^'A. (BI3) 
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Eq. (|B12[) is nothing but the continuity equation. 

The continuity equation derived above contains the 
current relaxation function. Thus, we need to next de- 
rive an equation of motion for this. The derivation goes 
along in a quite analogous way as that of the continu- 
ity equation. Specifically, to end up with an equation 
having ui4>j(q,ui), we will apply the following onto the 
Bethe-Salpeter equation, instead of SqG^ 1 ; 

5 j G' 1 {k;q,u)) ee i[5 G _1 (fc;g,a;), q»%] + . 

Since ui and q being sufficiently small, we will keep only 
its leading-order contributions; 

/ BF" \ 
^G- 1 ^;?,^)- (u ! -^-+2iFo)q„%-qi. (B14) 

Apply this onto eq. (|B1[) and take the sum over k, k' and 
the band indices. By way of this, we obtain the following 
constitutive equation; 

Qpl „ 

^-E^^W* 

k,k± 

x E^Xnft fc i; 9, w) E $ *7.rn ( fc i . 9, w)(B15) 



Aj and 7^(fc; q, w) are defined as follows; 



(B16) 



{<5G(fc;q,^) ■ G R '-\k+, (,+) ■ ^% ■ G R (k + ,fi+) + 
G A (k_,v-) ■ q„% ■ G A >-\k_,v-) ■ SG(k; g, w)JB17) 

Contrary to the continuity equation, this equation of 
motion contains the convolution between the 2PIR ver- 
tex function and the response function explicitly. This 
convolution part describes the interactions between the 
current relaxation function and the other types of relax- 
ation functions. We will linearize this convolution part 
with respect to relaxation functions in the following three 
paragraphs. 

To do this, note first the completeness of the 7- 
matrices, 

E hvM%W ( Bl8 ) 

^={o,l,.,S,15,— ,42} 

and that of the spherical harmonic function Yi m (Cl), 
00 1 

f(x(l)^J2 E YUfy^YTrnfiM*®), (B19) 
1—0 m— — l 0/ 




<tt« 




2 FIR 



m m 



FIG. 16: The convolution between the 2PIR vertex func- 
tion and the response function is replaced by the direct 
product between the relaxation kernels and relaxation func- 
tions, where we used the complete set for a function of hi, 
i.e. f(ki) = J2 k2 5 ( kl - fc 2)/(fe) = E fe2 E m u ™( fc i) x 
u m{k2)f(k2). Namely, u m constitutes the bare vertex part 
described by ">"-mark in the figure, while u* m constitutes that 
described by "<"-mark 



where Q denotes the normalized vector and En ' ' ' stands 
for the 2-dimensional integral over the angle-direction; 
J2q = 47r. Using these two completeness relations, 
we can decouple the convolution part in eq. (|B15[) into 
the sum over the countable numbers of modes (see also 
Fig. [HI); 



i^-Q^ + 2iF )0j - #o(«,w) 



k,k± 
15 oc I 

X EE E Y Uk)[%] Sl ^i m J\h\-,q,u). (B20) 

(J,=0 l—Q m— — l 

Namely, the fci-dependence of the response function is 
decomposed into the dependence on its radial coordinate 
I k\ I and that on the angle coordinate k\. At a price for 
this, the right hand side contains the summation over the 
azimuthal and magnetic quantum numbers, / and m. For 
each I, m and //, <j) lrn fl (x; q, u) is defined as follows; 

EE E [%]p a Kn$)*wp(xk,k>iq,u)J(p21) 

% k ' <x>P>1 

Observing this definition, notice that the x-dependence 
of (f>i m ,a(x; q, u>) and its to, g-dependence can be further 
factorized for small u> and q; 

l^lra.^l^) = 9lm,li{x)4>lm,li(q,u). (B22) 

This is because, for such small u> and q, the response func- 
tion appearing in eq. (|B21[) is dominated by the diffuson, 
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which depends only on lu and q\ 

x [f d ( g ,o;)] Qi5ii7i/3i Gf i5 (A : ;^ + )G^ i (A : '_, / i_).(B23) 



By taking the integrals over k and k' in eq. (|B21[) and 
keeping only the leading order in small to and q, one 
can actually verify this factorization for any I, m and (i 
(consult also the appendix. D for several examples.) 

Without loss of generality, we can assume that gi m A x ) 
thus obtained is normalized with respect to the integral 
over the radial direction; 



x 2 dxgi m ^(x) 



1. 



Then, corresponding 0/ mj/J (q, w) given in eq. (|B22[) be- 
comes the non-zero azimuthal number (7 ^ 0) generaliza- 
tions of the relaxation functions defined in eqs. (|B3IB4|) : 



Jim 



m(9,^) = E E [%] 0a Y *rn(k)^ (k,k'-,q,^2A) 

k,k' a,/3,7 



Thus, substitute eq. (|B22|l back into eq. (|B20| . Then, 
we finally obtain the constitutive equation, which is fully 
linearized with respect to these relaxation functions; 

(w-^p + 2iF )<pj(q,u) - q<p (q,uj) = Aj - - 

oo I 15 



(|fc'|)y im (fc')} (t>imAQ^)- (B25) 



/=0 m— — l ^,—0 k,k 



The second member of its right hand side described 
the 'mean-field' induced by other relaxation functions, 
into which the 2PIR vertex function is encoded. This 
situation is quite analogous to how the 1-particle irre- 
ducible function (i.e. self-energy) describes the interac- 
tion among 1-point Green functions. 

Due to this interaction, however, the constitutive equa- 
tion above also contains relaxation functions (jji m ^(q,(v) 
assigned to the higher order harmonic (I > 1) sector. 
Thus, to make the final coupled EOM's to be closed, we 
must also derive the EOM for all of these functions. This 
is, however, limitless. To make it tractable, we thus need 
to truncate interactions among these too many modes. In 
this paper, we will consider the interactions only within 
the "s-wave" sector. Namely, we will restrict the sum- 
mation over I, m and \i in eq. (|B25[) to the I = sector; 



,r, (?,<*>) = E E [%] /3a $ "7,7/3( fc I fc '; < ?,w). 
k,k f a,f3,*y 



These 16 modes in the s-wave sector further reduce into 
the 8 modes, when the rotational symmetry is taken into 
account. Namely, by noting that the response function is 



invariant under the simultaneous rotation in the pseudo- 
spin space and in the momentum space; 

a , a [KA Y ^*s„f3(k, k; q, u) J ss , [ulj 0fj , 

(J^Ee'W^', (B26) 
we can derive the following identity; 

15 , 45 j 

4;t 



15 - 45j 

E [%] aptao Ali = "TpE Q/3 ^(g,w)(B27) 



a=0 



Ai=0 

with: 

<t>a(Q,^)=^2 E [^a{q)] 0a ^aj,jl3{k,k';q,Uj), 
k,k' a,/3,7 

(■00, #5, #4, £45) = (70,75,74,745) 

(B28) 

By use of this equality, eq. (|B25|) turns out to consist 



label 


density 


T 


1 


label 


"current" 


T 




"0" 


7o 


+ 


+ 


"j" 


<7m 7m 






"5" 


75 


+ 


+ 


"5j" 




+ 




"4" 


74 






"4j" 


9m7m4 




+ 


"45" 


745 


+ 




"45j" 


?M ~ e A tI/ p7'^P 




+ 



TABLE III: Symmetry of eight modes in the s-wave sector and 
their symmetries under the spatial inversion 1 and the time- 
reversal T. Since 70, 74, 75 and 745 behave as a scalar quantity 
under the rotation defined in eq. (|23|l . we regard them as the 
"density" associated with the sublattice and spin degrees of 
freedom. Corresponding to these four types of density, we 
have 4 types of "current", which in turn behave as a vector 
quantity under the rotation. 

only of those 8 functions defined in eq. (|B28jl ; 

u-Q^ti - #0 + MjA<l,u)<t>a = Aj,(B29) 

" a=0,5,---,45j 



M a , b (q,v) = 2iFg5 ab + 



k,k' 



x [ 7a L (fc; q,A] pa U^(k, k'; q,u) [ 7b R (fc'; q,A] s ^0) 
7„' R (fc;g,w) above are given as follows; 

^(k;q,uj) = - X 

{SG(k;q,A ■ G«'- 1 (fc +! /i+) • v a (q) • G R (k+, t i+) + 

G A {k-, //_) ■ « (g) ■ G A ' _1 (&_,//_) ■ 6G(k;q,uj)}BZl) 
^(k;q,A = Va(q)9ooA\k\)- (B32) 
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M a .b(q, w) defined in eq. (|B30p generally appears in 
the EOM for the i; a -type relaxation function and plays 
role of the 'mean-field' induced by the Wb-type relaxation 
functions. Namely, this 8x8 matrix is nothing but the 
'self-energy' in the matrix-formed EOM's for the s-wave 
sector (see eq. (|B35|) ). Thus, we will refer to M(q,uj) 
as the relaxation kernel henceforth. Before deriving the 
remaining 6 constitutive equations, let us remark on the 
general property of this kernel. Observing eq. (|B30|) , no- 
tice that each element of this matrix becomes pure imag- 
inary, when its two arguments taken to be zero; 



M a ,b{q, w)* = -M a ,b(-q, -w). 



(B33) 



This can be directly seen from; 



{SG(k;q,u;)}* = -{6G(k;-q,-u J )} t . 
{U^{k,k'-q,u)Y = U$J a (k,k';-q,-u J ). (B34) 



The EOM's for the other 6 s-wave modes can be de- 
rived in parallel with that for the current relaxation func- 
tion. Specifically, we will begin with the Bethe-Salpeter 
equation applied by the following, instead of eq. (|B14|) ; 

8 a G~ l {k;q,Uj) = -[8oG~ l {k;q,w), V a (q)] + , 

with v a (q) taken to be 75,74, ^9^^%% respec- 
tively. Going through the same procedure as described 
so far, we will reach the constitutive equations for these 
remaining 6 modes. Combined with eq. (|B12[) and 
eq. (1B29[) . such equations consist of the following 8 by 
8 matrix-formed EOM's; 

[^(g,u)+%u)].feu) = 4u). (B35) 

<p(q,ui) and A(q,ui) have the eight components; 



£ 2 = 



Mi 



Mo 



with 



UJC 

ujd~ x + 2ie~ 1 — q 




2ie- 



LOC 



Mi 



(B38) 



M 5 ,o M 5 ,5 
M iA M iA5 

M45A M i5 A5 

MAZjAj M 4 bj,45j 



d ~^W- e " Sft ' <B39) 



By solving eq. (|B35|) . one can obtain the asymptotic 
expressions for the relaxation functions for small lo and 



00 (9, 



An 



r, <k(q,u) 



A Bo 



(B40) 



U) + iDq 2 ' ' oj + iDq 2 

where \Aq\ stands for the density of state (see eq. (|B13[) V 
The (renormalized) diffusion constant D used above and 
other coupling constants are expressed only in terms of 
the relaxation kernels estimated at u>,q = 0; 

M 5 j,5j 



D = i 



Bn = - 



M,,M:,, :,, - MijjMjfij J \ q , u =o 
2iF" + M 5 n 



M 



0.5 



(B41) 



(B42) 



q^LJ—O 



Eq. (|B41[) and eq. (|B30|) become the essential bulding 
blocks of our gap equation (see below) . 



2. gap equation and its solution 



W = [00, <f>j, 05, 05j, 04, 04j, 045 , 045 j] , 

A* = [A Q , Aj, A 5 , A 5j , A 4 , Ay, A 45 , Aisj], 



latter of which is defined as follows; 

A a {q,uj) = J-J2 T 4t(k;q,iu)]. 



K.(q,uj) and M(q,u>) are defined as follows; 



K = 



K.i 



/Ci 

6 fc 2 



M = 



Mi 
6 Mi 



ui —q 
—q ioc~ x 
ual- 1 + 2ie~ 1 



(B36) 



(B37) 



When the non-crossing approximation is employed for 
the 1-point Green function, the diffusion pole in eq. (|B40|) 
should be attributed to the ladder-type diagram whose 
long-wavelength expressions were already obtained in 
sec. IV. Especially, we have observed that in the sec- 
tion. IVB that the charge diffusion mode and parity diffu- 
sion mode equally dominates the diffuson in the massless 
case (Dl~ 2 > r t ~ po ), while the parity mode becomes inef- 
fective in the presence of the relatively large topological 
mass {Dl~ 2 < r t ~ po ) (seeeqs. (|106I98I102|I ). Correspond- 
ing to these two limiting cases, we will derive two types 
of gap equations and their solutions in this section. 



a. for m = case 

Let us begin with the zero topological mass case first. 
In this case, we will sum up eq. (|107p and eq. (|108p . since 
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h = fi- With use of eq. and 02,3,4 = 0, such a 
summand takes on a following form; 



f d ( q ,oj) lF ^ 



ao/4 



1 + Ti + 5i + 5 S 



(B43) 



In section. IVB, we have observed that the overall factor, 
ao/4, has the diffusion pole as in eq. (|98p . where its bare 
expression were calculated explicitly. Namely, by keeping 
track of the small q effect in eq. (|78p , we obtained the bare 
diffusion constant as in eq. . Instead of such bare ex- 
pressions, however, we will describe henceforth this ao/4 
in terms of the renormalized diffusion constant defined by 
eq. (|B41I) . Namely, we want ao/4 to be given by the relax- 
ation kernels, only to obtain the self- consistent equation 
for the diffusion constant. 

To do this, notice that relaxation functions for small 
q and uj are dominated by the diffuson as in eq. (|B23[) . 
Thus, by substituting eq. (|B43)) into eqs. (|B23IB3|) , we 
will first express the density relaxation function in terms 
of a / 4 ; 



(f>o(q,oj) ~ -647rio;A 2 ao/4 



(B44) 



The factor A 2 in the right hand side stems from the mo- 
mentum integral over k and k' in eq. (|B3[) : A is the ul- 
traviolet cutoff of the momentum-integral. Then, we will 
equate this with cf>o(q,uj) obtained in the step-(i), i.e. 
eq. (|B40p . By doing this, ao/4 is expressed in terms of 
relaxation kernels; 



64"7riaao/4 



1 



An 



A 2 uj + iDq 2 ' 



(B45) 



Namely, the diffusion constant D in the right hand 
side was already given by the relaxation kernels as in 
eq. (lB4il . 

Substituting this back into eq. (|B43|) . we obtain the 
asymptotic form of the diffusion; 



f d (a,w) 
ao/4 



4 

1 1 



1 + Ti + Si + S 2 
A 



The diffusion constant D in eq. (|B47p is now given by 
the relaxation kernels, via eq. (|B41|) . These relaxation 
kernels are in turn defined by the 2PIR vertex function, 
via eq. (|B30|) . The 2PIR vertex function is usually dom- 
inated by the Cooperon given by eq. (|B47|1 . at around 
k + k' ~ 0. As such, we will replace (approximate) the 
2PIR vertex function in eq. (|B30|) by this asymptotic 
form of the Cooperon, i.e. eq. (|B47|) . By way of this, we 
obtain closed coupled equations for the (renormalized) 
diffusion constant D: 



D 



M 



3J Jvl 5j,5j - Mr,,.jM 



14)1 



2 12 tt 2 L>A 2 



I L~ 1 <\k+k'\<l- 1 

a L (fc)uii - £ - si + s 2 } . 



' (57 



\k + k' 



(B48) 
d 3 kd 3 k' 

(B49) 



Since eq. (|B47|) is valid only for small |A; + k'\, we have 
imposed the additional constraint |fc+A;'| < l~ x into these 
integral variables. One might regard this upper limit as 
the mean-free path. We have already taken in eq. (|B49|) 
both uj and q to be zero. Thus, 7 L ' R (fc) in the right hand 
side stands for 7 L ' R (fc; g, uj) estimated there; 

t(k) = 7«(*;0 S 0), 
la( k ) = v a x goo, a (\k\), 



where v a was already defined in eq. (|B28|) . The normal- 
ized real- valued function goo. a {x) used above is given only 
in terms of Fq. For example, goo,j(x) is given as follows, 



.9oo j (x) 



4 1 



8 x 2 F'< 



3|F 2 



■},(B50) 



with its normalization factor N. 



3i 



2 8 nai A 2 uj + iDq 2 



When its hole line time-reversed, the corresponding 
Cooperon at small uj and k + k' is also derived; 

U coop (k + k\uj) 



aao/4 
4 

1 1 



(l + fx + Si + S 2 J ■ (B46) (see the appendix. D for its derivation). Thus, eqs. (|B48l 

IB49P constitute closed coupled equations for the diffusion 
constant. 

To solve this gap equation, notice first that the cou- 
pling between the current and the 75-type current is dis- 
connected in the massless case; Mj,bj = 0. This can be 
seen directly from 



1 - Ti - Si + S 2 
A 



- ¥ (l-f 1 -5i + 5 2 )(B47) 



1 - Ti - Si + S 2 ) q u [% 5 ] = 0. (B52) 



2 s m A 2 uj + iD{k + k'f 
where we used the following identities; 

[l®Sy] 7 , 7 Tl,a5,/3 7 [l«)S y ] /3/3; = -Ti >aSn ,p,, 
[1 ® Sy] 7 ^Si, a s,p 1 [i ® Syjpp, = -Sl <aSt yp>, 



which leads to D = i/M.j,j- As a result of this, eqs. (|B48I 
IB49P becomes linear in D; 



D 



2 7 ir 2 D A 2 



d 3 kd 3 k' 



L- 1 <\k+k'\<l- 1 



x gooAM) KF ') 2 -(^) 2 + fc 2 )-2(fc-g) 2 

\k + k'\ 2 ((F>) 2 - (Fq ) 2 - fc 2 ) 2 + 4(F>) 2 (F^') 2[ } 
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FIG. 17: A schematic plot of 2Dt 1 as a function of the 
chemical potential in the massless case. 



Using eqs. (|B50IIB5T|) . we can readily evaluate the mo- 
mentum integral in the right hand side of eq. (|B53[) . To 
do this, introduce a new integral variable q' = k + k' , 
so that dkdk' = dkdq' . Moreover, we approximate 
9oo,j(\k-q'\) in the integrand by goo,j(\k\), since goo,j( x ) 
is a slowly varying function in the scale of l^ 1 . These 
treatments give us the following expression for 2Dt~ 1 ; 



2Dt~ x = 1 + - 



2^ 



A 



-Ti 2 ' 

3^ 



(B54) 



where we used A = — 16Fq'A and F = fl+iT^ 1 . Observ- 
ing this expression, notice that the second member of the 
r.h.s. is nothing but the quantum correction to the dif- 
fusion constant, which basically corresponds to the AWL 
correction to the conductivity. 



b. for m^O case 



In the presence of relatively large topological mass, i.e. 



DV 



< T. 



topo : 



the parity diffusion mode becomes the 



high-energy degree of freedom. As such, only the first 
term in eq. (|106|) contributes the diffuson. With eq. (pj)) . 
such diffuson is given as follows; 

f-WW = -£^{l±^(i + f 4 )(i + s x ) 



i-t 2 



(1 + T 4 ) (fx + S 2 ) + t(f 2 + T 3 ) (1 + Si) },(B55) 



where s and t are defined by ao4, a x and 023; 



a 4 



l + r 



-3ai 



1-t 2 



-a 23 = st. 



Contrary to the previous subsection, the tensor-part of 
the diffuson depends on the model-parameters through 
the "tensor-form factor" t. As such, we will employ 
in this case not only the diffuson constant D but also 
this tensor-form factor t as the "mean-field parameters" , 
which should be self-consistently determined. In other 
words, both of them should be given by the relaxation 
kernels, as in eq. (|B41[) . 



To do this, we will first calculate both the density re- 
laxation function 0o and the sublattice density relaxation 
function </> 5 , by the use of eq. (1B55[) . Namely, we will sub- 
stitute eq. (|B55I) into eqs. (|B23IB28|) . only to obtain these 
two functions in terms of sfa and t first. The relaxation 
functions thus calculated read as follows; 



<t>o(q,u) 



-32naA 2 isfa, 



<p5(q,ui) ~ 32iraA 2 istfa. 



(B56) 
(B57) 



Then, we will equate eqs. (IB56IB57I) with the first two 
members of eq. (|B40I) respectively. By way of this, fa 
and t can be given in terms of the relaxation kernels; 



t 



1 1 



\Ao\ 



32na A 2 ico - Dq 2 ' 



= -B = 



2iif 



(B58) 
(B59) 



By substituting these two back into eq. (|B55j) . we can 
express the diffuson only in terms of the relaxation ker- 
nels. When its hole-line time reversed, the corresponding 
Cooperon is readily derived; 

U coop (k + k',cu) ~ 

, . n ^° ;n2 {(i + t 2 )(i + T4)(i-^) 

2 9 7ri A 2 lj + iD(k + k') 2 L y /v ' 

-(1 - t 2 )(l - fk) (Ti - S 2 ) + 2t(f 2 + fa) (i - Si) ]{B60) 

The tensor-form factor t and the diffusion constant D ap- 
pearing in the right hand side above are already given by 
the relaxation kernels, via eq. (|B59p and eq. (|B41[) . Such 
relaxation kernels are given by the 2PIR vertex function 
(see eq. (|B30I) ) . Thus, as in the previous subsection, we 
will approximate the 2PIR vertex function by eq. (|B60p . 
In terms of this substitution, we arrive at a closed coupled 
equation for the diffusion constant D and the tensor-form 
factor t, whose explicit expressions are given in the ap- 
pendix. C. 

When solving this gap equation, we can see how the 
quantum correction to the diffusion constant behaves as 
a function of m and fi. Several limiting values are sum- 
marized in Fig. 1181 Especially, in the zero mass limit i.e. 
m = 0+, the solution of the gap equation reduces to a 
following simple function of Fq = /I + ir~ 



-1. 



lim 2£>t~ 1 = 1 



1 l- x -L- 



±U 2 



12 



A 



f=2, (B61) 
3 /i 



Comparing this with eq. (|B54[) . one can easily see that the 
quantum correction to the diffusion constant is actually 
half of that for m = case. 

The discrepancy between eq. (|B54I) and eq. (|B6ip is 
responsible for the Cooperon term associated with the 
parity diffusion mode, i.e. the second member of the 
right hand side of eq. (|112|) . To see this explicitly, note 
first that a 2 , 03 and 04 appearing in eqs. (|113H114"|) reduce 
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12 A 



FIG. 18: A summary of several limiting values of the quan- 
tum correction of the diffusion constant in the presence of the 
topological mass. Eq. (|B61[) corresponds to the green shaded 
region. 

to zero in the limit of m — > 0+. Then, eqs. (|113I114[) in 
this limit read as follows; 

lim U^ = -\(i-f 1 -S 1 +S 2 ) 

+ f 4 -(l + f 1 -S 1 -S 2 )}, (B62) 
lim U^ = -\(l-f 1 -S 1 + S 2 ) 

m->0+ L v ' 

_r 4 .(i + fi -&-&)}. (B63) 

Observing these expressions, notice that the second mem- 
bers of both eq. (|B62p and eq. (|B63[) are totally ineffective 
in the current- type relaxation kernels; 

{f 4 . (i + fx - & - S 2 )} aSnp q,[%] Sy = 0,(B64) 

I 



{f 4 ■ (1 + 2\ - Si - S 2 )} aSi0 q, [% l5 ] h = Q(B65) 

The consequence of these two equations are two-fold. 
The second equation in combination with eq. (|B52[) leads 
Aij.5j = first. Thus, we have D = j^— again, which 
indicates that Mjj originated from eq. (|B62|) and that 
from eq. (|B63[) contribute to the quantum correction in 
an additive way. Eq. (|B64[) moreover indicates that these 
two quantum corrections have the same magnitude and 
sign. In other words, the quantum correction derived in 
Appendix. B2a, i.e. eq. (|B54|1 . can be divided into two 
parts; 

2-D|m=o T 1 — 

1 I- 1 - L- 1 t- 2 + \T? 1 r 1 - L- 1 t- 2 + \ji 2 
+ 12 A r- 2 + ±7I 2 + 12 A T- 2 + \jl 2 ' 

Each of these two quantum corrections is originated from 
eq. (|B62[) and eq. (|B63[) respectively. Since we have al- 
ready ignored eq. (|B63[) for the m ^ case, the resulting 
solution has only single 1/12, as in eq. (|B61[) . 
APPENDIX C: MEAN-FIELD EQUATION FOR 
m CASE 

The mean-field equation for the diffusion constant D 
and t in the presence of the finite topological mass m is 
given as follows; 

D _ t _ 2iF£ + M 5fi 
MjjM 5 j,5j - M 5 j,jMj, 5 j ' M:,.r, 

with the relaxation kernels M. a .b being given by D and t 
self-consistently; 



x it (.*>)} Pa {(i + ^ 2 )(i + r 4 ) (i - Sx) - (1 - < 2 )(i - t 4 ) (fx - s 2 ) + 2t(f 2 + f 3 ) (i - &) } Q(5i7/3 {7a R (fc')} 57 - 



1 

Note that 7 L ' R (fc) above are previously defined; functions (Fq and F5) and the tensor-form factor t; 

ti(k)=-{5G(k-,0,0).G R '-\k,»).v a .G R (k,ri x 2 + {(\F \ 2 + \F 5 \ 2 )-t(F a *F 5 +F Q FZ)} 

+ G A (k, f i)-v a -G A ^ 1 (k,n) -6G(k;0,0)}, 5oo,o(^) oc |(a + »6) 2 - a: 2 | 2 

7a R (fc) ^f) a x ff00 , a (|fc|), to 2 - {t(|F | 2 + |F 5 | 2 ) - (F*F 5 + F F 5 *)} 

\{a + ib) 2 -x 2 \ 2 

with w a for a = 0, 5, j and 5j given in eq. ([B28]> . g 00 ,a(» , v f ^0 ~ ^5 , 8 afr(i^ - f^)a; 2 \ 

used in 7 R (fc) are given in terms of the 1-point Green 9oo,j{x) oc j ^ + ^ 2 _ ^ + g ^ + ^2 _ x2|4 j>, 
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5oo, 5j (x) oc 
with (a + ib) 2 



8 abjFg - tF({)x 2 

\{a + ib) 2 -x 2 \ 2 ~ 3 \(a + ib) 2 -x 2 \ 4 



of eq. (|D1|) as follows; 



Fl 



APPENDIX D: DERIVATION OF goo,a( x ) 

Starting from the Bethe-Salpeter (BS) equation for the 
response function, we have derived in the section Bl the 
EOM's for the various types of relaxation functions. Such 
coupled EOM's have two features; they are closed and 
linearized with respect to the relaxation functions. Be- 
cause of these two features, we can solve them for the 
relaxation functions. Out of this solution, we can relate 
the renormalized diffusion constant with the 2PIR (two- 
particle irreducible) vertex function. This relation in turn 
becomes an essential building-block of the self-consistent 
loop of the diffusion constant (see appendix. B). 

To obtain such linearized EOM's, we need to reduce 
the convolution part between the 2PIR vertex function 
and the response function into the simple product be- 
tween relaxation kernels and relaxation functions. For 
this purpose, we have introduced the completeness in the 
space of the integral variable, say y or w, associated with 
this convolution; 



J2 Ua ^ =5(y-w). 



Namely, by use of this, any convolution in principle can 
be decomposed into a simple product; 



dyf{- ■ ■ ,y)g(y,- ■ •) 



dyf(- ■ ■ , y)u a {y) ■ / dwu* a (w)g(w, ■ ■ •) 



In the current context, f(--,y) corresponds to the 2PIR 
vertex function, while g(w,--) to the response function. 
Therefore, J dy ■ ■ ■ f(--,y)u a (y) corresponds to the relax- 
ation kernels, while J dwu*(w)g(w, ••)••• does to the re- 
laxation functions (see also Appendix. Bl). The trade-off 
for this decomposition is therefore the sum over infinite 
(but countable) numbers of modes specified by a. 

To be more specific, we did this decomposition sys- 
tematically, based on the completeness relation of the 
7- matrices and the spherical harmonic function Yi m (Cl); 

dwu* a {w)g{w, ■••)—> 

E E [%]pY: m m a ...,..A\k\k ■ • •). (di) 

k a -P 

As such, the momentum integral only over the angle- 
direction, k, is taken, while that over its radial direc- 
tion, |fc|, is not taken. As for the convolution with re- 
spect to this radial direction, we simply replace the In- 
dependence of <&,,,(\x\x,- • •) by some real-valued func- 
tion <7...(|x|). Namely, we rewrite the right hand side 



k a 3 

= <WW)E E [%] YC m {k)^ a ...,.. (k, ■ ■ -).(D2) 

k a,/3 

Let us justify this treatment of the radial direction. 
In the response function, $> a s n p(k,k';q,uj), q and lo are 
associated with the external momentum and frequency 
for the bosonic degrees of freedom. We can take these 
two to be small, as far as the relaxation functions for 
the long wave-length and low-energy region is concerned. 
Then, such a response function is usually dominated by 
the diffuson T d (q, lo) (see eq. (jllOp 1 ). As a result of this, 
the |fc|-dependence in the left hand side of eq. (|D2p and 
its q, w-dependence can be factorized at the leading order 
in small q and lo; 

EE E M^mW^TM^M';^) 

k' k a >Pn 

= 9lmA x )Yl E [%] /3 a Y *m(k)^a 7 ,y0(k,k'-, ,q,(j) 
k : k' a.^^ 

= 9im,ii(x) x 4>i mifl (q,u). (D3) 

To see this factorization more explicitly, one can take 
the following steps: (i) substitute the asymptotic tensor- 
form of the diffuson into eq. (|110[) and the left hand side 
of eq. (|D3[) , (ii) take the integral and the sum over k, 
k', a, P and 7 in eq (ID3|) , and (iii) retain the leading 
order in small q and lo. By way of this, one can reach 
the factorization given in the right hand side of eq. (|D3|) 
with a specific gi m ,n{x). 

For example, let us follow these prescriptions in the 
case of zero topological mass case. Observing eq. (|B43I) . 
notice first the following relation; 



r 



(D4) 



for /1 = 1,2,3. Using this, one can readily checked that 
the diffuson in this case turns out to be proportional to 
the unit matrix, when its right hand side is traced out; 

E ra,5, 7 /3(9, W )|f5=0 E ^£( fc +> V+)G£y{k'_, fl-) 

5,7 fc',e 



E 



\Fo\ 



>./2 



lu + iDq 2 ^ \F 2 - k 



/2|2 



As such, to obtain the normalized function goo. a (%) m 
the massless case, we have only to calculate the following 
quantity up to the leading order in small to and g; 



E 

k,k',S,a,/3 



1 



[v a {q)] 0o ^ a s,5fi{k,k'\q,Lo) 

E GUk+^+)Gt p {k-^-)[v a {qj\ 



iu - Dq 2 ^ ^ 
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For example, taking the current component as v a (q) 
above, we have; 

1 ( 1 8 F 2 k 2 1 mil \ 



Observing the right hand side, one can then convince 
oneself of eq. (|D3|1 . Moreover, the normalized real- valued 
function g 00 j(x) will be obtained as in eqs. (|B50IB51|) . 
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